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^ Pattern formation in biological, chemical and physical problems has re- 

ceived considerable attention, with much attention paid to dissipative sys- 
tems. For example, the Ginzburg Landau equation is a normal form that 
describes pattern formation due to the appearance of a single mode of insta- 
bility in a wide variety of dissipative problems. In a similar vein, a certain 
^ "single-wave model" arises in many physical contexts that share common 

pattern forming behavior. These systems have Hamiltonian structure, and 
Ch the single- wave model is a kind of Hamiltonian mean-field theory describ- 

2 ing the patterns that form in phase space. The single-wave model was 

originally derived in the context of nonlinear plasma theory, where it de- 
scribes the behavior near threshold and subsequent nonlinear evolution of 
unstable plasma waves. However, the single-wave model also arises in fluid 



> 

in 

\^ mechanics, specifically shear-fiow and vortex dynamics, galactic dynamics, 

the XY and Potts models of condensed matter physics, and other Hamil- 

O 

tonian theories characterized by mean field interaction. We demonstrate, 

m 

by a suitable asymptotic analysis, how the single- wave model emerges from 
a large class of nonlinear advection-transport theories. An essential ingre- 
i dient for the reduction is that the Hamiltonian system has a continuous 

. 5^ spectrum in the linear stability problem, arising not from an infinite spatial 

domain but from singular resonances along curves in phase space whereat 
^ wavespeeds match material speeds (wave-particle resonances in the plasma 

problem, or critical levels in fluid problems). The dynamics of the contin- 
uous spectrum is manifest as the phenomenon of Landau damping when 
the system is stable, and when the system becomes unstable, an embed- 
ded neutral mode appears within the continuum. The single-wave model 
captures the dynamics of the embedded mode if the system is perturbed so 
as to make the mode unstable, or if perturbations destroy that mode and 
the system remains stable. As a consequence, the model describes a range 
of universal phenomena: for a bifurcation to instability, the model fea- 



tures the so-called "trapping scaling" dictating the saturation amplitude, 
and the "cats-eye" structures that characterize the resulting phase-space 
patterns. In the stable situation, the model offers one of the simplest de- 
scriptions of Landau damping, and illustrates how this can be arrested by 
nonlinearity. Such dynamical phenomena have been rediscovered in differ- 
ent contexts, which is unsurprising in view of the normal-form character of 
the single-wave model. 
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I. INTRODUCTION 

A. Perspective 

The dynamics of forming patterns in biological, chemical, and physical systems has received considerable 
attention over recent years (see e.g. Cross and Hohenberg, 1993). In the simplest settings, one deals with 
spatial patterns that form as a result of the appearance of a small number of "modes" of instability, such 
as a convection cell in a fluid contained between differentially heated plates. These arc the eigcnmodes 
of the linear stability problem and provide a natural set of tools to describe the forming spatial patterns. 
In particular, by decomposing perturbations about the unpatterned equilibrium into the eigenmodes and 
exploiting the center-manifold reduction, the full system of governing equations can be approximated by a 
low-order set of amplitude equations near the onset of linear instability (e.g. Crawford, 1991). 

The center-manifold reduction works because, in spatially bounded, dissipative systems, the eigenvalues 
of the normal modes compose a set of distinct points on the spectral plane, and instability arises when one 
of these modes becomes unstable. The purpose of the present article is to summarize extensions of some of 
these ideas to non-dissipative systems that suffer instabilities that cannot be described in low-dimensional 
terms, and for which the patterns that form occur naturally in the associated phase space. These systems are 
Hamiltonian systems that have a neutrally stable continuous spectrum of eigenvalues in the linear stability 
problem. Moreover, this spectrum is not the result of the problem being couched in an infinite spatial 
domain, but due to singularities in the linear eigenvalue problem that arise due to a resonant interaction 
between wave-like perturbations and the background equilibrium state; the singular resonances occur along 
curves in phase space whereat wavespeeds match material speeds. Such systems are also commonplace, 
and include ideal plasmas, inviscid fluids, self-gravitating stellar systems, and various other Hamiltonian 
models characterized by mean-field interactions. The plasma problem is that of electrostatic instability in 
an ideal, single-species plasma (van Kampen, 1955), which was the original building block in understanding 
a multitude of more complicated instabilities in plasma theory and fusion science. The fluid problem is the 
instability of inviscid shear flow or vortices (Case, 1960), which are also key problems in fluid mechanics. 
Another mean-field Hamiltonian model of interest is that with the XY interaction (see Chaikin and Lubensky, 
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1995), which models spin-spin interaction and has been used to describe e.g. superfluid helium and hexatic 
liquid crystals (see the 'Hamiltonian Mean Field' model of Campa et al., 2009 for an extensive treatment). 

The mathematical issue which underlies our discussion is that, although instability is described by an 
exponentially growing discrete eigenmode, the mode detaches from the continuous spectrum when it becomes 
unstable. Consequently, at onset, the distinguished mode cannot be isolated from an infinite number of 
other, singular eigensolutions representing the continuum (the center-manifold reduction relies on such a 
spectral gap) . This aspect of the problem fundamentally affects the weakly nonlinear description, which, as 
a consequence, has to proceed down a relatively novel and much more tortuous pathway. 

This pathway leads to the single- wave model of O'Neil, Malmberg & Winfree (1971), following terminol- 
ogy in plasma physics, or simply the "single-wave model", for short. Tliis model is a Hamiltonian normal 
form that describes the transition to instability described above in a wide variety of physical systems of fluid, 
plasma, and other disciplines. The model captures the dynamics of the patterns that form in the phase space 
of the system due to the growth and nonlinear saturation of the unstable mode detaching from the contin- 
uous spectrum. All of the physical models for which the single-wave model is appropriate are hyperbolic 
Hamiltonian systems with Hamiltonian characteristic equations. Thus the patterns occur in a conventional 
phase space witli canonically conjugate coordinates, and consequently have a Liouville theorem and other 
Hamiltonian properties, all of which help to recognize the circumstances that lead to the single- wave model. 

For the physical systems we treat as examples, the patterns have two characteristic features: there is a 
relatively smooth, wave-like pattern ocupying most of the domain, coiipled to a more complicated, finely 
scaled structure appearing over a localized region surrounding the resonance of linear theory. The larger- 
scale wave pattern reflects the shape of the embedded neutral mode. The finely scaled structure is generated 
by the interaction of that mode with a locally resonant packet of the continuous spectrum, and takes the 
form of either localized "holes" in the plasma distribution function (see, e.g., Eliasson and Shukla, 2007), 
"cat's-eye" vortices in the fluid shear flow, and 'magnetization' zones in the XY model. Explicating the 
universal nature of these patterns as described by the single-wave model is the main goal of this article. 



B. Single-wave model preview 

The single- wave model is a Hamiltonian normal form that captures the dynamics of Landau damping in 
stable situations, and the nonlinear saturation of unstable modes. In light of its origin, the model must be 
of Hamiltonian and transport form, and be defined in a phase space that we take, for simplicity, to be two- 
dimensional. The basic dynamical variables are a density-like variable Q{x, y, t), and an amplitude variable 
A{t). The density Q depends on the coordinates {x,y), describing the phase space over the resonant, finely 
scaled region. The amplitude A{t) depends purely on time, and represents the amplitude of the embedded 
neutral mode; given that mode's eigenfunction, one can reconstruct the large-scale wave pattern outside 
the resonant region. In addition, the model contains parameters that must be adjusted to characterize the 
physical system of interest. A special case of the single- wave model is compactly given by following: 

Qt + [Q,S] = 0, (1) 

ip = ^e'^ + A*e-'^ , (2) 

iAi = (Qe--), (3) 

where S = — cp^ 

[f, a] = fx9y - fvdx , and {■) = — J dy dx ■ (4) 

and fx ■= df/dx, etc. Equation (1) is the transport equation for Q that resembles the Vlasov or Euler 
equation, with dependence on the variable ip through £, an energy for the characteristic 'particle orbit' 
equations. From the second of Eqs. (3) it is clear that <f represents a wave of amplitude A, with its single- 
wave spatial dependence upon the variable x being obvious, while equation (2) gives the temporal evolution 
of A, which is driven in an integral manner by the rearrangements of Q throughout the entire resonant 
region. For simplicity, the full parameter dependence of the single-wave model has been chosen in (l)-(3) 
to furnish the simplest possible form of the equations; thus no parameters appear explicitly. 
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FIG. 1 A solution of the single- wave model in its normal form. Parameter values are k = and 
7 = — 1 (cf. Sec. V). The system is initialized with A{0) = and Q{x,y,0) = —y and then 
kicked into action by adding a forcing term O.lt exp(— 200t'^) to the mode amplitude equation (see 
Balmforth and Piccolo, 2001 for further details of the computational scheme). Panel (a) shows 
three snapshots of Q{,x,y,t) as densities on the {x,y)-plane. Two wavelengths are shown for the 
final snapshot to emphasize the periodic geometry. Panel (b) shows the time evolution of the mode 
amplitude, A{t); the stars indicates the times of the snapshots in (a). 



The solution of (l)-(3) is presented in Fig. 1, which shows three snapshots of the density Q{x, y, t) on the 
(x,y)-phase plane (panel a), and the mode amplitude A{t) as a function of time (panel b). For the earliest 
snapshot, Q{x,y,t) has not been significantly rearranged, but at later times the Q distribution twists up 
into a cats-eye pattern. This is an example of the universal pattern seen in many physical systems described 
by the single-wave model, and is typical of vortex formation in fluid mechanics and electron hole formation 
in plasma physics. 

The system of Eq. (2) is a Hamiltonian field theory given in Poisson bracket form as 



Qt = {Q,H} and 



{A,H} 



where the Hamiltonian is given by 



H 



dy / dxSQ 
-oo Jo 



(5) 



(6) 



and the Poisson bracket { , } is 



{F,G} 



OF dG dF dG 



2n \ dA dA* dA* dA 



oo /'27r 

dy I dxQ 



5F SG 
6Q'6Q 



(7) 



functional derivatives such as 6F/6Q are defined in Sec. II. This Hamiltonian structure is inherited from 
that of the various 'parent' models for which the single-wave model is an encompassing normal form. 
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C. Overview 

In Sec. II we describe a large class of Hamiltonian 2+1 field theories that contain the basic ingredients 
for the phenomena described by the single-wave model. In this section, we also establish our notation and 
provide examples of specific physical systems, including shear flow from fluid mechanics, Vlasov theory from 
plasma physics, and the XY model from condensed matter physics. We use, mainly, these three examples 
to help make our general points throughout the paper. We note that the mean field model with the XY 
interaction has been called the 'Hamiltonian Mean Field model' by Campa et al. (2009) and collaborators. 
We refrain from this terminology since all of the models treated in this paper are in fact examples of 
Hamiltonian mean field models. Instead, we will refer to this example as the XY model. 

As mentioned above, an essential feature of single-wave dynamics lies in the nature of the linear theory 
of the parent models. In particular, to understand single-wave phenomena it is necessary to understand 
the CHH transition to instability that it describes. This transition has been developed and rediscovered 
independently in many fields, and since this linear dynamics is crucial to our development, we review this 
in an historical context in Sec. HI. We first describe the linear theory for our general class of Hamiltonian 
models and then consider the early history due to Rayleigh and Kelvin, followed by the developments in 
plasma physics, stellar dynamics, and in the XY model. A main purpose of producing a normal form theory 
is to unify discoveries made in disparate fields; the nature and considerable amount of rediscovery and 
interchange in several areas emphasizes the need for a review of the single-wave model. 

As is common in normal form theory, asymptotic calculations must be performed and the situation for 
the single-wave model is particularly intricate. Thiis, in Sec. IV we outline the nonlinear dynamics that one 
must understand to appreciate the single-wave model. In particular, we first demonstrate how conventional 
weakly nonlinear theory is unable to describe the transition to instability and fails to capture the scale of the 
nonlinear saturation level. That level is referred to as trapping scaling in plasma physics, and we provide 
a historical commentary on the origin of this scaling, as well as parallel development in fluid mechanics 
regarding what is termed critical layer theory. These descriptions motivate the asymptotic sequences used 
in the derivation of the single-wave model. 

With the intuition gained from Sees. HI and IV, the stage is set for the derivation of the singe-wave model, 
which we do in Sec. V. The derivation involves an asymptotic expansion featuring both the discrete embedded 
mode and the continuous spectrum. We provide a general derivation that produces the Hamiltonian normal 
form embodied in the single-wave model. We also observe that, in some specific physical problems, a certain 
degeneracy occurs that necessitates a degenerate form of the single-wave model, a form shared by top-hat 
profiles in the XY model and some waterbag distributions for the plasma and stellar dynamics problems. 
Thus we delineate a norm form and a degenerate form of the single-wave model. Next, we discuss the 
Hamiltonian nature and conservation laws of the single- wave model that are inherited from the Hamiltonian 
parent model. Finally, for later use, we outline the leading-order forms expected for different types of 
dissipation. 

The linearized version of the single-wave model is arguably the simplest system capturing the transition 
to instability of a discrete mode embedded in the continuous spectrum, as well as the disappearance of such 
a mode and Landau damping. In Sec. VI wc describe this linear theory in detail, including the effect of 
dissipation, so as to address issues of structural stability of our model. Various patterns emerging from the 
single-wave model are described in Sees. VII and VIII. In Sec. VII we discuss patterns of the normal form: 
first where the pattern arises out of instability and second where the pattern arises from nonlinc;arly forcing 
a stable situation. In Sec. VIII we provide a similar analysis for the degenerate form, considering in detail 
the XY model. In Sec. IX several variants of the single-wave theme are discussed. Finally, in Sec. X, we 
conclude and discuss the universality and limitations of the single-wave model. 



II. A CLASS OF HAMILTONIAN THEORIES AND SPECIFIC EXAMPLES 

We now describe a general class of 2-|-l Hamiltonian mean-field theories and outline several specific 
examples that we discuss in more detail in later sections. 
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A. A class of Hamiltonian field theories 

The class of 2+1 Hamiltonian field theories (Morrison, 2003) possesses a single independent variable, 

l^{q,p,t), which is a dcnsity-likc variable that depends on the independent phase space variables z (^jp), 
coordinates for some phase space Z, as well as time. We suppose the density satisfies a transport equation 
of the following form: 

f + [C,f]=0, (8) 

where the 'particle' Poisson bracket is defined by the usual expression [f,g] = fqQp — Qqfp, where fq := df/dq 
etc., and the quantity £ is an energy- like quantity that we call the particle energy. If (8) were a Liouville 
equation, then £ would be a given function of z, and we would have a linear theory. However, we are 
concerned with mean field theories where £ depends on in a global sense. Thus, our systems are governed by 
nonlinear (quasilinear) partial integrodiffcrcntial equations. Such equations arise, for example, by truncation 
of BBGKY-like hierarchies, which results in particular functional dependencies of the particle energy on the 
density. We do not pursue such an approach, but postulate that the particle energy arises from a total 
energy (Hamiltonian) functional of the form H[C,] = Hi + H2 + H3 + . . . , where the one-point energy. Hi, 
the two-point, H2, etc. are given as follows: 



Hi[C]= Jd'zhi{z)C{z,t).^ 



H2[c] = l 1^'^ I d'z'az,t)h2{z,z')az',t), (9) 



2 



z Jz 



with the generalizations to iJs, H^, etc. being obvious. The quantities hi and ft,2, the interaction kernels, are 
left unspecified for now. But, we do suppose the two-point interaction possesses the symmetry hiiz^z') = 
hiiz! , 2;). The particle energy is obtained from the field energy by functional differentiation: 



^ := ^ = /ii + ^d^z' h2{z, z') C(/, t) , (10) 



bH_ 



where the functional derivative is defined as usual by 6H = (fiz 5( SH/6(. It is not difficult to show that 
H[Q is a constant of motion under the dynamics of (8). 

Equation (8) with £ = 5H/5Q is a Hamiltonian field theory (Morrison, 2003). However, since there is 
only one field variable Q the theory is not of a canonical form, but possesses a noncanonical Lie-Poisson 
bracket description (for review see, for example, Morrison, 1998 and Marsden and Ratiu, 1999) given by 



{F,G} = jd^zC 



5F_ 5G_ 



(11) 



Noncanonical brackets of this form, which depend explicitly upon the variable C unlike conventional Poisson 
brackets that only depend on (functional) derivatives of the canonical variables, express the Hamiltonian 
form of matter in terms of Eulerian variables. Such noncanonical Poisson brackets were introduced in the 
context of magnetohydrodynamics in Morrison and Greene (1980), and the specific form of (11) was given for 
the Vlasov-Poisson system in Morrison (1980) and the two-dimensional Euler equation in Morrison (1982). 
Using (11) the equations of motion are obtained in the form 



5H 



-[C,^]. (12) 



Associated with the preceding Hamiltonian form are various constants of motion. The Hamiltonian H[Q 
itself is one such constant. In addition, there are the Casimir invariants, 



C[C]= Jd'zCiO, (13) 



where C(C) is an arbitrary function, which arise from degeneracies in the Poisson bracket. 
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Other important invariants arc momenta, P[C], that arc Hamiltonian dependent. Momentum invariants 
generally arise from translation symmetries that in the present context would be determined by the form 

of the interaction kernels /ii, /12, This is how the strong version of Newton's third law is built into the 

n-body problem. This idea can be generalized in various ways, hut we do so by observing that our system 
conserves momentum if there exists a canonical transformation of the phase space Z 

z = {q,p)^z:= {6,1) 

such that in the new particle coordinates z := {0,1), the interactions hi, /12, etc. have upon composition 
with z{z) one of the following two forms: 

hioz = hi{I), h2o{z,z') = h2{I,l',\e-9'\) (14) 



or 

hioz = 0, h2o{z,z')=h2{\I-I'\,\e-9'\). (15) 

In the first case 

p[c] = jji^ziaz). (16) 

is conserved, while in the second case we have two kinds of translation invariance and thus two components 
of the momentum 

Pi[Q= [ cPziaz) and P2[C] = [ d^zO C{z) . (17) 
Jz Jz 

All examples treated in this review will be of one of these types and three-point and higher interactions will 
not be considered. 

Associated with each mean-field model is a Hamiltonian n-body problem. This follows by assuming a 
Klimontovich type of distribution 

n 

C = Y.5{z-z\t)), (18) 

i=l 

with z^ = {q^,Pi), i = 1, 2, . . . , n, and substituting (18) into (8) and seeking a weak solution. Alternatively, 
the n-body problem arises upon effecting the functional chain rule on = f{z^,z'^, . . . ,z"), where / is 
defined upon inserting of (18) into F, giving 



d SF 
5^ If 



to establish a mapping between the canonical and noncanonical Poisson brackets. Both methods yield the 
Hamiltonian n-body problem 

z' = [z\H], (20) 

where self-energy terms are removed from H. The existence of this Hamiltonian n-body problem reinforces 
why our patterns occur in a phase space in the mean-field theory, rather than a configuration space. One 
consequence of this is that the characteristic equations of the single-wave model possess Liouville's theorem 
on conservation of phase space volume (here area). 

All of the examples treated here arise from equilibria that only depend on I. For this reason we set 
C,{9, 1, t) = R{I) + p{6, 1, t) and then when a choice of R is made, p{9, 1, t) represents the main dynamical 
variable. We further assume that 9 is an angle-like variable, so that the phase space is periodic in ^ e [0, 27r], 
and I gV with V equal to [—1, 1] or (— cxd, cxd), depending on the example. Upon substitution of ( = R + p 
into £, both of the forms of (14) and (15) can be written as follows: 



£[R + p\=£[R]+£[p] =: h{I) + ^{e,I), 



(21) 
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with 



^{e, I) = icp:= f di' f d0' K{i, I', \6 - 0'\) p{e', r, t) , 

Jv Jo 



(22) 



where h and K are determined by hi and /i2- Thus the governing equations become 

Pt + [R,^] + [p,h + ^]=0 (23) 



or 



Pt + npe-R'^e + [p,'^] = 0, (24) 

where [f,g] = fegi — gefi and r2(/) = ft'. The quantity //(/) corresponds to a Hamiltonian for an integrable 
system and the coordinates {9, 1) are action-angle coordinates for the characteristics of this system. 

Equation (24) will serve as our starting point in subsequent analyses of our examples. In addition for all 
the examples we consider, the interaction is assumed to have zero mean, i.e. 

/ d0' K{I,I',0-e') = O. (25) 
Jo 

Of note is the case where K is independent of I and /'; such a K will be referred to as a spatial kernel. 
B. Notation 

By tradition, different independent and dependent variables have been used in different contexts. For 
example, in plasma physics, the phase space coordinates, representing position and velocity, are often written 
as (x, v). In fluid mechanics, both coordinates denote position, just as the example of Sec. LB used (.x, y). To 
be a little more definite, we will always describe the first independent variable as coordinate-like or angular, 
while the second will be momentum-like or an action, as for the general independent variables, {q,p) and 
{0,1), defined above. Moreover, throughout, the following two shorthands will be used: 

[/, g] = fq9p- fq 9p (26) 
{f) = ^jjpj^\f{q,p), (27) 

where {q,p) denotes any 'conjugate' pair and, as above, subscript denotes partial differentiation. 

When there is a proliferation of subscripts, we will use e.g. dx to denote d/dx etc. In all cases our 
coordinate-like variable will be periodic and eventually scaled to be 27r-periodic, while the momentum 
domain V is example specific. When V = (— oo, oo) we will denote it by K. 



C. Examples 

1. Shear flow 

Our example of shear flow considers the evolution of disturbances to a channel flow or to perturbations 
of differentially rotating vortices, as illustrated in Fig. 2. Thus, we begin with Euler's equation in two- 
dimensions, in which case the field ( corresponds to the scalar vorticity field. In cylindrical geometry the 
coordinates for Z would be (0,1) = (0,r^/2), the usual polar coordinates; however, we are principally 
interested in channel flow with a periodic boundary condition and so we adopt {9,1) = {x,y), where the 
spatial coordinate x is equivalent, after a simple scaling, to the "angle" 0. Thus, for shear flow, the action- 
like coordinate / is equivalent to y, the cross stream coordinate, the equilibrium channel flow speed in the 
x-direction is U{y), and /C is the inverse of the two-dimensional Laplacian for this geometry. We then set 
C = R + p = U'{y) + uj{x,y), which here physically corresponds to divorcing the background equilibrium 
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FIG. 2 Sketches of the geometry of our inviscid shear flows and ideal plasma. In the channel flow 
in (a), the flow proceeds in the x-direction, and shears in the y-direction according the equilibrium 
profile, U{y). The unstable disturbances considered in this article take the form of weakly unstable 
modes that have a pronounced effect within a localized "critical layer" which surrounds a line, the 
"critical level" along which a neutrally stable wave propagates at the same speed as the mean flow. 
In (b), the circular variant of the problem is shown: a differentially rotating vortex with equilibrium 
rotation rate, rio(r). The critical ring, or co-rotation radius, now locates the critical region. The 
sketch for the plasma problem in (c) displays the total distribution function, F{v) + f(x,v,t), as 
a surface over the {x, t>)-plane; the location of the wave-particle resonance and the surrounding 
critical region are indicated. For the XY model example, shown in (d), the distribution function, 
f{6,I,t), develops a cat's eye pattern around the resonance at I = 0. 



flow from an evolving perturbation described by a relative vorticity field uj(x,y,t). Thus, (24) for this case 
is 

uJt + Ulo, - f/'V, + [t/^, cj] = , (28) 

where (22) is equivalent to 

1^ = '4'xx+ ^yy , (29) 

and 'ip{x,y,t), the streamfunction, is identified with $. Equations (28)-(29) are dimensionless, the channel 
width and maximum equilibrium flow speed acting as the relevant units. The flow is bounded by impermeable 
walls at y = ±1, and so ipi^i il; t) = 0> ^ind periodic in < a; < L, where the domain length is L. 

Thus, to summarize, the shear-flow problem is recovered from the system of Sec. II. A by the following 
identifications: 

0^x, I^y, nil) o Uiy), R{I) ^ U'iy), 
p LO, and (formally) K ^ A^^ . 
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2. The Vlasov-Poisson plasma 

In the Vlasov-Poisson system, the plasma is described by an electron distribution function (phase space 
density) depending on a spatial coordinate, ,r, and speed, w, with a stationary neutralizing ion background. 
The relevant equilibrium distribution function is spatially homogeneous and determined by a profile, F{v). 
Disturbances about this equilibrium are described by the perturbation distribution function, f{x,v,t), and 
electrostatic potential, ip{x,t), that satisfy the Vlasov equation, 

ft + vU + 'PxFv + /] = 0, (30) 

and Poisson equation, 

V>xx = dv f{x, V, t)- Nb, (31) 

where Nb is the neutralizing backgroimd (i.e. Nb = if))- Here, the electron mass, charge, and an assumed 
thermal velocity scale are used to nondimensionalize the equations. Thus, lengths are measured in units of 
Debye length and time in units of the electron plasma frequency. We assume periodic boundary conditions 
in X, and we demand that / vanish as w — > ±oo. 

Thus, the plasma problem is recovered from the system of Sec. II. A by the following identifications: 

61 o a;, I ^v, n{I) o V, R{I) ^ F{v), f, 
ifi -H- and (formally) K ^ = \x- x'\ . 

3. XY model 

Our XY model captures essential features of long-range spin-spin interaction and has been studied in a 
variety of contexts (Campa et ai, 2009; Chaikin and Lubensky, 1995). The phase space coordinates (6,1) 
and basic dynamical variable f{9,I,t) are equivalent to those for the Vlasov-Poisson system of Sec. II. C. 2, 
and the model differs only by a change in the interaction: 

ft + Ife - ^eRi - <^efi = 0, (32) 

$ = -- /d/ / do' f{e',i,t) cos{e-e'). (33) 

^ Jm J a 

Thus, the XY model is the system of Sec. II. A with the following identifications: 

e^e, 7 0(7) o 7, R{I) o R{I), 

f, and K ^ -tt"^ cos(6l - 6') . (34) 



4. Other models 

There are many models with physical content that fit into the same general Hamiltonian framework 
described here. One example is the Jeans equation for stellar dynamics (obtained by changing the sign 
of the interaction in the Vlasov equation and removing the zero net charge condition). Other examples 
include various quasi-geostrophic models for describing ocean and atmospheric dynamics (Pedlosky, 1987). 
In shear flow, the vorticity defect model of Balmforth et al. (1997) (obtained by simplifying the relationship 
between the vorticity and the streamfunction) is also applicable. Interesting related examples where the 
underlying characteristics correspond to integrable n-body problems are the Cologero-Moser system (Illner, 
2000; Moser, 1975) and a model introduced by Smereka (1998). 
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III. LINEAR THEORY AND HISTORICAL PERSPECTIVE 

The linear theory of fluid and plasma systems with continuous spectra has a long history. Here we 

present a general theory in Sec. III. A for the bifurcation to instability, followed by discussions of the early 
developments in fluid mechanics in Sec. III.B, plasma physics in Sec. III.C, and other systems in Sec. III.D. 



A. General eigenspectra; bifurcation to instability through an embedded neutral mode 

We begin by substituting p = p{I, t)e''^^ (with wavenumber k) into Eqs. (22) and (24) of the general 
formulation of Sec. II. A. After linearizing the latter, we arrive at 

Pt = -ikflp + ikR'KkP, tkP = [ dl' kk{I, I')p{I', t) , (35) 

Jt> 



and 



'0 

Special cases of Kh of interest are the following: 
shear flow: 



Kk{I,I')= / deK{e,I,I')e"'r (36) 
Jo 



Kk = \ 

Vlasov plasma: 



C sinh[fc(J+l)]sinh[A-(j'-l)] -r _ 1 < / < // 
' ksinh[2k] — — 



\ sinh[fc(J-l)]sinh[fc(/' + l)] .r 7' < 7 < 1 



Kk = -k' , 

XY - interaction : 

Kk = -Sk,i - 6k,-i , 

with 5i,j being the Kronecker delta symbol. 

The first term on the right of (35) corresponds to free-streaming motion of fluid elements or particles 
along each level of /. In the Fourier representation, the free-streaming dynamics is analogous to a continuum 
of uncoupled oscillators with a range of frequencies 0(/) indexed by /. The shearing action of the basic flow 
corresponds to phase mixing due to the frequency spread. This can be seen particularly clearly for cases in 
which the second term on the right of (35) is absent, leading to 

p = p(7,t)e"=^=p(/,0)e'^[«-^^(^)*l 
= p{1, 0) j dce"=(^-^*) 5{n{I) - c) , (37) 

where C denotes the range of values of f2(/) that correspond to X>. This solution models the tilting over of 

the initial perturbation by differential advcction; the lines of constant phase are given hy 9 — Q{I)t = const. 
The final piece of (37) further emphasizes how the solution cannot be written in a form that is separable in I 
and t unless one formulates it as an integral superposition of singular modes of the continuous spectrum with 
delta- functions for eigcnfunctions (Dirac, 1927; Eliasscn et ai, 1953; van Kampen, 1955). In other words, 
one can identify the free-streaming operator in the term, Q{I)p$, as creating the continuous spectrum. 

The second, interaction term in (35) acts as a perturbation of the free-streaming operator. However, 
with certain restrictions on the form of /C/j(/, /'), this perturbation of the operator docs not change the 
existence of the continuous spectrum (for rigorous treatments see Degond, 1986; Hagstrom and Morrison, 
2011; Kato, 1966; Morrison, 2003), and only introduces the possibility of additional discrete eigenmodes. 
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This was shown by cxphcit construction of the solution of the initial-value problem by Case (1960) and 
Rosencrans and Sattinger (1966). 

To find the discrete normal modes introduced by the interaction term, we set p{I,t) = p(/)e~*'^'^*, 
furnishing 



R'r] 

P= ] 

and 



P^'qZTc' '^^^kP (38) 



If c is complex, this is a regular integral equation. But in the problems of interest, the integral on the 
right-hand side of (39) becomes singular when c is real, which is symptomatic of the continuous spectrum. 
Moreover, for c = Cr + ici and Cj ^ 0^ , 



ri{I) 



where means the singularity of the integral is evaluated by taking the Cauchy principal value, and the 
Tj 's denote all values of I for which r2(/) = Cr (Gakhov, 1990). Wc may satisfy the real and imaginary parts 
of (40) if we choose R'{Ij) = for all the /j's, and take rj{I) to be real. This choice is straightforward if 
f2(/) is a monotonic function, as we assume later when we derive the single- wave model. Then, Ij — ^ J* 
and R'{I*) = 0(/*) — = 0. Non-monotonic equilibrium with certain symmetries can also satisfy this 
condition (such as the jet portrayed below, in Fig. 6, which is symmetrical about the midline). In either 
circumstance, the special mode that is embedded within the continuum and which is also the limit of the 
complex eigenmodes as Cj 0^ is the smooth solution of (40). The distinguished modes are therefore 
associated with the extrema of the equilibrium profile, R{I) (inflection points of the velocity proflle, ft = 
U{y), in the fluid problem). 

Note that, if the kernel is independent of the action coordinates, what we call a spatial kernel, then 
/C — > /C(6') and the JC^'s and therefore the eigenfunctions r] are constants. Hence, (39) reduces to an explicit 
dispersion relation, 

D(c,fc) = l-/C, / d/-^^=0. (41) 

For the plasma problem, D[c, k) is the plasma dispersion function, and its relatively simple form allows one 
to extract a number of general criteria for determining the transition to instability, results such as those of 
Penrose (1960) (see Sec. III.C). Similar results have recently been outlined for the XY model (Hamiltonian 
mean field model) by Chavanis and Delfini (2009) as well as vorticity defects by Balmforth et al. (1997), 
although they amount to a minor generalization of those provided by Penrose, which are standard in plasma 
physics (Krall and Trivelpiece, 1973). 



B. Rayleigh's problem 

In the late 1800s, Rayleigh explored a variety of fluid stability problems motivated by experiments with 
smoke jets and flames. As vividly illustrated by cigarettes and chimneys, smoke jets often rise initially 
uniformly, but then begin tortuous undulations and meanders. Rayleigh's theoretical explanation for the 
breakdown of the unidirectional shear flow in the jet amounted to a linear stability analysis of (28) and (29). 
For exponentially growing, wave-like, normal-mode solutions with ip oc e''^^^"'^*), the system reduces to the 
so-called "Rayleigh equation," 

(42) 
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This equation follows from (39) on recalling that the kernel is Green's function for the fluid problem, and 
then recasting this equation in its differential form. 

Rayleigh found explicit examples of unstable modes for inviscid shear flows in which U{y) took the 
form of a sequence of broken lines (corresponding to a stacked set of layers with constant mean vorticity). 
For such profiles, U"{y) = everywhere except at the line breakages, which allows for an explicit closed 
form solution of (42). His results for these broken-line profiles motivated Rayleigh to continue further to 
look for a more useful stability criterion for general flow profiles, and he formulated what is now called the 
"inflection-point theorem" : we multiply (28) by tp* and integrate across the channel. After a little algebra, 
involving integration by parts and separation of real and imaginary parts, we arrive at 



Thus either Cj = (and there is no exponentially growing mode), or the integral must vanish, which demands 
that U" change sign somewhere within the channel. In other words, a sufficient condition for stability is 
that U {y) have an inflection point. Further results that extended Rayleigh's theorem were presented later by 
Balmforth and Morrison (1999); Fj0rtoft (1950); Howard (1964); and Rosenbluth and Simon (1964). These 
criteria can be argued to be equivalent to energy stability criteria (Balmforth and Morrison, 2001), akin to 
Dirichlet's theorem of classical mechanics (see, e.g, Morrison, 1998). 

Rayleigh's equation (42) has singular points for real wavespeeds where U{y) = Cr- These singularities 
typically prohibit the construction of neutral waves (e.g. Balmforth and Morrison, 1999), and occur at the 
so-called "critical levels" for which neutral waves propagate at the same speed as the background flow. These 
locations play a special role in the single-wave theory described later, and are examples of the resonances 
mentioned in Sec. I. Note that the critical-level singularities appear only in the normal-mode problem: as 
is clear from our discussion in Sec. III.A, the linear initial-value problem itself is completely regular, and it 
is only by forcing solutions into the normal-mode form that one produces the singularity. 

The singularities in Rayleigh's equation, and their prohibition of neutral waves led Kelvin to object quite 
strongly to Rayleigh's theory: "This disturbing singularity vitiates the seeming proof of stability contained 
in Lord Rayleigh's equations." Kelvin's view was that Rayleigh's singularities reflected a breakdown of 
linear inviscid theory. At the critical levels, "the motion has a startlingly peculiar character" , nonlinearity 
could not be ignored, and any slight disturbance created "a cat's-eye pattern of elliptic whirls". In fact, 
Kelvin was convinced that viscous shear flows were linearly stable, as illustrated by his analysis of Couette 
flow {U{y) = y). Kelvin was also probably prejudiced by Reynolds' experiments in pipes, which had been 
undertaken at roughly the same time, and suggested that there was a finite threshold in the amplitude of 
perturbations required to generate turbulent motion. Kelvin argued further that inviscid flows were strongly 
unstable in the sense that the slightest perturbation would generate cat's eye patterns along any critical 
levels ("let one or both bounding-surfaces be infinitesimally dimpled in any place and left free to become 
plane again... Hence the interior disturbance essentially involves elliptic whirls. Thus we see that the given 
steady laminar motion is thoroughly unstable, being ready to break up into eddies in every place, on the 
occasion of the slightest shock or bump on either plastic plane boundary" ) . Kelvin evidently felt that there 
were some fundamental problems with the inviscid formulation. 

Despite Kelvin's criticisms, Rayleigh remained convinced that inviscid fluid theory was a useful guide, 
and that shear fiows could be unstable, perhaps guided by his experiments with fiames and smoke jets. 
Rayleigh contended Kelvin's objections by pointing out that the inviscid linear analysis was perfectly valid 
for exponentially growing modes - in this case c is complex and there are no critical level singularities. 
However, he did admit that his proof of stability for flows without inflection points only concerned expo- 
nentially growing instabilities and did not exclude algebraically growing instability ("Perhaps I went too 
far in asserting that the motion was thoroughly stable; but it is to be observed that if c be complex, there 
is no disturbing singularity. The argument, therefore, does not fail, regarded as one for excluding complex 
values."). In fact, it was not until the 1960s that algebraic instability was to a large degree ruled out 
for shear flow when researchers followed Landau's Laplace transform treatment of the plasma problem (cf. 
Sec. III.C) and solved the linear problem as an initial-value problem (Case, 1960; Dikki, 1960; Engevik, 
1966; Rosencrans and Sattinger, 1966). 

Rayleigh also criticized Kelvin's assumption that viscous flows were stable: "Lord Kelvin arrives at the 
conclusion that the flow ... is fully stable for infinitesimal disturbances.... Naturally it is with diffidence 
that I hesitate to follow so great an authority, but I must confess that the argument does not appear to me 




(43) 
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demonstrative." Nevertheless, some of the inviscid profiles that Rayleigh aimed to establish were unstable to 
explain unsteady motions in experiments were actually shown to be stable by the inflection-point theorem. 
Rayleigh's impression was that the theoretical problem lay in taking the inviscid limit ("it is possible that 
the investigation in which viscosity is altogether ignored is inapplicable to the limiting case of a viscous 
fluid when the viscosity is supposed infinitely small"). Indeed, as it is now known, the viscous version of 
the linear stability theory yields the celebrated Orr-Sommerfeld equation, rather than (42), which can have 
unstable modal solutions even for flow profiles that are stable in inviscid theory (Landahl, 1986). 



C. Plasma oscillations and Landau damping 

In the middle of the 20th century the analog of the disturbing singularity discussed by Rayleigh and 
Kelvin's arose in plasma physics in the context of linear analysis of the; Vlasov equation; the singularity now 
appears at the wave-particle resonance where the particle velocity matches the speed of a neutral plasma 
wave (oscillation). To observe this directly, let / = /(i>)e''^^^~'^*' -|-c.c., and then linearize equations (30) 
and (31) to obtain 

which clearly presents the singularity at the wave-particle resonance, = c for c = real. Unlike for 
Rayleigh's equation the Vlasov problem has a spatial kernel and, thus, we can proceed much further ana- 
lytically and derive a dispersion relation for the wavespeed c: eliminating f{v) yields 

(cf. Eq. (41)). Complex conjugate solutions for c can be straightforwardly found from this dispersion relation 
for certain equilibrium distribution functions, F(v)] again, the mode with > is unstable. Even without 
specific choices of F{v), Nyquist methods can be brought to bear on the dispersion relation (45) to construct 
explicit stability criteria (Penrose, 1960). 

The singularity of plasma oscillations is seen in a different way in (45): this dispersion relation contains 
a singular integral with a branch cut along the real axis of the complex c plane, which again identifies the 
continuous spectrum (for rigorous spectral theory see Degond, 1986, and Hagstrom and Morrison, 2011). 
On using the Plemelj relation (Gakhov, 1990) with c ^ c, -|- iO^, we find that 

1 f F'(v) 

F'{c^) =0 and ^-^27^^ — = ° • (^^) 

In other words, as the equilibriiim profile F(v) is adjusted so that the complex conjugate pairs move towards 
the real axis of the complex c-plane, they limit to a special neutral mode with c = c*. By virtue of the first 
condition in (46), the special neutral mode is a smooth solution (the singularity of l/(w — c^) is cancelled 
by the zero of F'{v) in the numerator of the eigenfunction) , even though it is embedded in the continuous 
spectrum. The embedded mode can be interpreted to lie at a location within the continuous spectrTim 
where the modal "signature" changes sign, consistent with a generalization of Krein's theorem (Hagstrom 
and Morrison, 2011; Morrison, 2000) of finite-dimensional Hamiltonian dynamics (Krein, 1950). 

For a family of equilibria indexed by a control parameter a, the embedded mode determines the stability 
boundary on the (fc, a)— plane because it is the limit of an unstable complex eigensolution. By way of an 
example, consider the family of equilibria given by 

F{v)=e-''' +ae"^''-''°^\ (47) 

where a and Vq are additional parameters. This parameterized family consists of a primary Maxwellian 
with a superposed "bump-on-tail" . Via Penrose's criterion one can establish that there are no distinguished 
neutral modes associated with the main peak of the distribution near v = 0; instability can only arise when 
the bump on the tail generates further extrema. If the bump amplitude, a, is the main control parameter, 
the stability boundary is denoted by the curves a = a{k), as illustrated in Fig. 3. The minimum of this 
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FIG. 3 (a) Stability boundary on the {k, a) plane for parameters values of (47), a = 4 and vq = 2. 
The star and cross mark the values of a for which the profile is plotted in (b). The star denotes 
the marginally stable profile for = 1 (o, ~ 0.182); the cross marks an unstable profile with a 
slightly larger bump (a = 0.23). The inset in (b) shows a magnification of the bump-on-tail region, 
with the shaded area illustrating the critical region surrounding the minimum of the marginally 
stable profile over which the distribution function is wound up by the unstable mode. This region 
is narrower than the velocity spread of the bump. 

boundary (a « 0.103 for k w 0.62, a = 4 and vq = 2) occurs when an inflection point first appears in F{v) 
in the vicinity of the bump, and thereafter splits into two extrema with a further increase of a. However, 
if the domain is periodic in x, as considered here, there is a minimum value of k determined by the longest 
spatial wavelength, and so the transition to instability can take place elsewhere on the stability boundary, 
as illustrated in Fig. 3. 

For stable equilibria below the stability boundary, the preceding arguments rule out all regular eigen- 
modes, and in particular neutral plasma oscillations. This physically unsatisfying conclusion led Landau 
(1946) to attack the linear initial- value probk^ni using Laplace transforms in order to uncover unambiguously 
how waves evolved once introduced into the plasma. A direct consequence of Landau's solution was the 
prediction of the phenomenon that is now called Landau damping: the exponential decay of the electric field 
(or potential ip{x,t)) that arises due to the phase mixing of particles. 

The exponentially decaying contribution to the electic field is often refered to as a quasi-modc (or Landau 
pole), as it is not a true normal mode of the linear eigenproblem. To uncover these objects, and the general 
long-time behavior of the solution, one inverts the Laplace transform solution, which takes the form of a 
complex integral over the usual Bromwich contour {e.g. Krall and Trivelpiece, 1973). The Landau poles 
can be picked up as the leading-order long-time contributions by deforming the integration contour through 
the branch cut of the continuous spectrum. One thereby emerges on a diff'erent Riemann sheet of the non- 
analytical dispersion hmction D{c, k) (which appears in the denominator of the Laplace transform solution), 
allowing one to encircle and pick up residue contributions from singularities due to zeros of the dispersion 
relation on a non-physical Riemann sheet. These "fake eigenvalues" are the quasi-modes. 

Despite Landau's reputation as one of the eminent theoretical physicists of his day, the relatively math- 
ematical nature of Landau's predictions was not convincing to many plasma physicists. Several questions 
arose. For example. Landau damping required that the energy contained in the electric field must be re- 
moved, which appeared to some to violate energy conservation (but see Morrison, 1994, 2000; Morrison and 
Pfirsch, 1992 where it is shown explicitly how energy is conserved and energy is transferred). Others argued 
that the damping effect was a mathematical artifact arising from inconsistencies in the formulation (Allis, 
1959; Allis et al, 1963; Ecker and Hoelling, 1963; Weitzner, 1963). On the other hand, in three celebrated 
papers Bohm and Gross (1949a,b, 1950) argued for the existence of Landau damping outside of Vlasov the- 
ory. The issue was not settled until experimental verification of Landau damping was reported by Malmberg 
and Wharton (1964) (sec also Malmberg et al, 1966), nearly twenty years after Landau's original paper. 
These authors demonstrated experimentally that the damping vanished when the resonant particles were 
removed. They did this by removing the tail of the distribution for electrons streaming in one direction but 
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not the other, and observed Landau damping only for the waves propagating in the direction with resonant 
particles. 

Prom the mathematical perspective the story of Landau damping has continued. Steps toward an early 
rigorous justification of Landau damping by Backus (1960), who considered the finear limit, were recently 
extended to the nonlinear regime by Mouhot and Villani (2011). Their proof requires a condition on the 
initial amplitude in order for the damping to behave essentially like the linear theory. However, as we will 
see later, the linear results fail for finite amplitude distiubances (see Fig. 10 of Sec. VII. B). The nonlinear 
arrest of linear Landau damping remains an open problem from a formal mathematical perspective. 



D. Violent relaxation, XY model, and phase transitions 

Closely connected to Vlasov-Poisson dynamics is the relaxation of stellar systems, which can be described 
by the Jeans equation (Binney and Tremaine, 2008). which is identical to the Vlasov-Poisson system but for 
an attracting interaction (the sign of /C is switched) and the removal of the charge neutrality condition. In 
a famous paper Lynden-Bell (1967) prompted a discussion of this problem by arguing that such relaxation 
proceeded in two phases: in the first, the ensemble of stars behaves like a continuum, relaxing rapidly from an 
arbitrary initial condition towards some quasi-steady state (the "violent relaxation"). That state consisted 
of a structure, or pattern, in phase space characterized by increasingly wound-up filaments. Eventually, 
the ever decreasing scales precipitate the failure of the continuum approximation, and herald the onset of a 
second phase in which the discreteness of the system plays an essential role. 

More recently, Lynden Bell's violent relaxation was revisited in the context of n-body Hamiltonian 
dynamics (e.g. Campa et al., 2009) with the XY interaction. Like gravity, this interaction is attractive, 
and the stellar and XY models has many common features. In addition to the consideration of far-from- 
equilibrium initial conditions, recent explorations with the XY model also studied equilibrium distributions 
with low-amplitude perturbations. In particular, stability theories were presented to identify equilibria with 
linear instabilities (Antoni and Ruff'o, 1995; Campa et al, 2009; Inagaki and Konishi, 1993; Yamaguchi 
et al., 2004), and numerical simulations explored the resulting phase-space patterns. Specific attention was 
focussed on the nature of the phase transition occuring when the linear instability appeared, and whether 
statistical techniques could be exploited to predict integral measures of the emerging patterns (such as the 
magnetization; see Antoniazzi et al. (2007) ; de Buyl et al. (2009) ; Campa et al. (2009) ; and Yamaguchi et al. 
(2004)). Unsurprisingly, the phenomenon of Landau damping was also rediscovered (Campa et al., 2009). 

The search for instabilities in the XY model closely follows the Vlasov-Poisson problem: after introducing 
the normal-mode form, f{9,I,t) = /(7)e'^^~'^*^ -|- c.c, into (32) and (33) and linearizing, we arrive at the 
explicit dispersion relation, 



which is the current version of (41) since this problem has a spatial kernel. All Fourier modes with k ^ ±1 
are stable. 

Two simple examples for the equilibrium distribution, i?(/), that allow further analytical progress are 
provided by the "top-hat" profile, 




(48) 



i?(/) = |[e(/ + i)-e(7-i)], 



(49) 



where 6(a;) represents a Heaviside step function, and the Gaussian, 




a 



e 2 



(50) 



(cf. Campa et al., 2009; Yamaguchi et al., 2004). For the top-hat profile, R' consists of a pair of delta- 
functions and so the integral in (48) is performed immediately to furnish 
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Thus, when a > 1, there is an unstable mode. This example is analogous to the "waterbag" model of the 
plasma problem (see e.g. Berk and Roberts, 1967). Moreover, both are equivalent to Rayleigh's broken- 
line profiles. For the Gaussian, the dispersion relation (48) is closely analogous to the plasma dispersion 
relation with a Maxwellian equilibrium distribution function (Fried and Conte, 1961). The neutral stability 
condition can be computed on observing that marginally stability occurs for c = 0, in which case a = ttc = 1. 
Moreover, once again, for a > 1, there is an unstable mode. 

Note that the instability of both distributions in (32) and (50) for a > 1 is quite different from Vlasov- 
Poisson theory, for which both equilibrium profiles would be stable. The key difference between the dispersion 
functions in (45) and (48) is the reversal of the sign of the integral, which arises because the XY interaction 
is attractive. The sign difference permits the distributions to become unstable in the XY model, a situation 
that is forbidden in Vlasov-Poisson theory because of the positivity of the distribution function. Likewise, 
the gravitational attraction of the stellar dynamics problem also reverses the sign of the integral in (45); the 
instability of a singly peaked distribution that results is a form of the classical Jeans instability, as explored 
by Fujiwara (1981). 



IV. NONLINEAR DYNAMICS: SCALING, SATURATION AND CRITICAL LAYERS 

Before proceeding to our derivation of the single-wave in Sec. V, we here motivate the necessity and 
features of the derivation. In particular, we show why conventional weakly nonlinear analysis fails, and 
review both the plasma physics and fluid mechanics literature. 



A. Conventional weakly nonlinear analysis 

The power of bifurcation theory lies in its ability to construct nonlinear states analytically. Convention- 
ally, the center manifold reduction (e.g. Crawford, 1991) is often the method of choice for dissipative systems, 
whereas Birkhoff normal form theory is applied to Hamiltonian systems (e.g., Bryuno, 1988; van der Meer, 
1985). 



1. Generalities 

Conventional weakly nonlinear theories aim for an amplitude equation that captures the dynamics of the 
normal mode that bifurcates to instability at a critical threshold in a control parameter. The derivation of 
the equation can be cast in the form of a regular perturbation expansion, in which one opens with a neutrally 
stable normal mode balanced on the stability boundary, and then kicks that mode into an unstable action 
by suitably adjusting the control. The relevant amplitude equation is dictated by the symmetries of the 
system in question; for the problems of interest here, the equation is a certain Hamiltonian normal form. 
This bifurcation is characterized by asymptotic scalings for the mode amplitude, the perturbation to the 
control parameter, and a redefinition of time. 

More specifically, and in the general scheme of the full model given at the end of Sec. II. A, the idea 
is that the equilibrium profile, R{I), is one of a family parameterized by a constant parameter, a. After 
a suitable choice of domain size and scaling, we focus on the situation in which the normal mode with a 
wavenumber fc = 1 is unstable if a > Oq, and all the other normal modes are stable. By switching the 
angular frame of reference we can zero the frequency of that mode (i.e. by transforming coordinates to a 
frame rotating at the mode's frequency). We then set a = ao + £^a2, where e <§; 1 is the small parameter that 
we use to organize the perturbation expansion. The modification to a embodied in e^a2 creates a change 
in the equilibrium profile that destabilizes the; othcTwise neutral mode. Correspondingly, the relatively slow 
growth of the mode is captured by scaling time. Thus, we set 

d d 

T = et, Tr,^e—, p{e,I,t) ^ p{e,I,T) etc. (52) 
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The remaining asymptotic scalings are encoded in the following sequences: 

i?(7) =i?o(/)+£^i?2 (/) + ... , (53) 
p = sp,{9,I,T)+e^P2{e,I,T) + ... , (54) 
^ = ICp = eM^,I,T)+s^ij2{^,I,T) + ... , (55) 

where Vm = ^Pm for 'm = 1,2, The expansion of p opens at order s, even though the equilibrium is 

modified only at O(e^), which corresponds to the "Hopf scaling" in the terminology of Crawford (see, e.g., 
Crawford, 1995). 

If one substitutes the sequences above into the system (22) and (24), and gathers terms of equal order 
in e, one arrives at a hierarchy of equations to solve sequentially for (prmV-'m)- The leading-order equations 
determining {pi,tpi) turn out to be the relations satisfied by the distinguished neutral mode. That solution 
has an undetermined amplitude, A{T). To fix this quantity, we proceed to the next orders. At 0(£^) the 
pair {p2,'>p2) capture the nonlinear generation of the first harmonic of the neutral mode (with wavenumber 
k = 2) and a correction to the angular average (with wavenumber k = 0), but A[T) yet remains unknown. 
Finally, at 0{e^), the requirement that the pair (p3,V'3) be bounded and periodic in 6 demands that A{T) 
satisfy a solvability condition, which amounts to the desired amplitude equation. The construction is quite 
standard, although for the current problems it contains a key flaw. We bring out the nature of this flaw and 
illustrate the expansion in a slightly simpler setting in Sec. IV.A.2. 



2. Weakly nonlinear expansion for spatial kernels 

By way of illustration, we perform a weakly nonlinear analysis of the special case of our general model 
for spatial kernels for which the kernel is independent of the actions, I and /'. As remarked earlier, explicit 
dispersion relations are available for such kernels, and both the Vlasov and XY models furnish examples. 

For further brevity, we consider the specific situation in which the domain is symmetrical in /, the 
equilibrium profile R{I) is an even function, and 0(7) := h'{I) is an odd function. One can then explore 
the situation in which an instability arises due to the appearance of a (zero-frequency) normal mode with 
the form of a standing wave, eliminating any need to transform into the rotating frame of the mode. After 
the rescaling of time, T = et, the governing equations become 

epT + ^P0 - ^sRi - '^^epi = (56) 

$= /" de' [ dl' K{0-0')p{0\I',T), (57) 
Jo Jt> 

where recall V denotes the range of /. It is also convenient to take the angular average, 
of the first of these equations: 

epT = (Wpe)i. (59) 

Next, we pose the asymptotic sequences, 

R = Rq + E^R2 , 

p = epi -I- e^p2 H , 

^ = 6^1+ £^$2 + • • • , (60) 

(cf. (53)-(55)). Introducing these into (56) leads, at order e to 

Pi = ^^, $i=A(T)e'^+c.c., (61) 
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after adopting the marginally unstable mode (with wavenumber fc = 1) as the leading-order solution. The 
relation in (57) then gives the marginal stability condition, 



1-Krldl^=0 (62) 
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(cf. the dispersion relation (41) with c = 0). 
At order e^, we find 



Hence, on introducing (61), 



\ Q J J 

Substitution of (64) into (57) yields 



i}p2e = Roi^ie - Pit - ^lePii- (63) 



1 (Roi\ 



(^V^+cc). (64) 



where 



h it). ■ 

Here we have used Kq = Q and that the symmetry of Ro{I) and ri(/) implies dl R'^/^"^ = 0. In principle, 
we should also add a second term to this solution with the form, j42(T)e'^ + c.c, which denotes a correction 
to the k = 1 neutral mode, but is not needed to determine the amplitude equation for A. 
At 0(£^) we also arrive at the first non-trivial result from (59): 

P2T = {^WP2 + ^2ePi)i = (^^) i\A\^)T. (67) 

To integrate this relation, we must specify an initial condition. If A{0) = Aq and /2(/, 0) = 0, corresponding 
to initializing the system with the equilibrium plus a low-amplitude normal mode, then 

— f -^0^ 

Finally, we arrive at the 0{e^) equation 



r2=\^-^)^{\A\'-\Aon (68) 



Pae = ^{Roi^ae + R2i^ie - P2t + ^2epu + ^iep2i)- (69) 

Isolating the contribution of the right-hand side to the first Fourier mode, e'^, and then applying the second 
relation in (56) leads to the amplitude equation, 



dI^ = JA + T\A\^A, (70) 



with 



and 



2{Ki - K2) Jv ^ 



Roi \ 1^ / Rqi 



(72) 
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Note that the integrals that appear in the coefficients of the amphtuclc equation do not exist unless a number 
of the derivatives of Ro{I) and R2{I) vanish at / = 0; this is reflective of critical- level singularities in the 
weakly nonlinear expansion. In fact, as shown by Crawford and Jayaraman (1999), the singularities appear 
at all orders of the expansion, becoming progressively more severe as one descends through the asymptotic 
hierarchy. The only instance in which one can justify the asymptotics is if all derivatives of Rq and i?2 vanish 
at the wave-particle resonance / = 0, which places severe limitations on the equilibrium profile. Notable 
examples of such situations, however, include Rayleigh's broken-line shear-flow equilibria, the waterbag 
models of plasma physics, and the top-hat profile of the XY model. 

3. Top-hat profile of the XY model 

For the XY model with the top-hat equilibrium profile, f2 := 7, D := K, Ki = —1, K2 = and 

Rji= "^[6(1+1) -6(1-1)], (73) 

with j = and 2. If we also assume that |^o| 1, the integrals can be evaluated to furnish 

ATT = a2A+l\A\''A. (74) 

Note that the sign of the second term ensures that nonlinearity is destabilizing in this example. That 
is, the bifurcation is subcritical, and an unstable steady solution branch given by \A\ = ^—la-i exists 
for a < ao = 1 (where the unpatterned equilibrium state is stable). A similar type of bifurcation arises for 
broken-line shear-flow profiles in the vorticity defect model of (Balmforth et al., 2012) (where Kk := (2fc)~^). 

B. Trapping scaling versus Hopf scaling in the Vlasov problem 

In the 1960s and 70s, several attempts were made to continue the analysis of linear electrostatic instability 
into the nonlinear regime. Some of the key results were obtained by Frieman et al. (1963); O'Neil et al. 
(1971); Onishchenko et al. (1971); and Simon and Rosenbluth (1976). One of the striking features of these 
articles is a disagreement regarding the level of saturation of the instability: the standard way to estimate 
the saturation level is in terms of the amplitude; of the electric field disturbance, Agat, relative to the distance 
to the stability boundary, as measured by a suitable control parameter, e (such as e = a — ttc, where Oc 
denotes a point on the stability boundary, for the example of Sec. III.C). Two characteristic scalings 
were hypothesized: "trapping scaling" (O'Neil et al, 1971; Onishchenko et al, 1971), for which Agat ^ e^, 
and Hopf scaling (Frieman et al., 1963; Simon and Rosenbluth, 1976) with Agat ~ e^^^) which is the usual 
saturation level for strongly dissipative instabilities and equivalent to that used above in Sec. IV.A. 

Trapping scaling was the essential ingredient in the "single- wave model" of Drummond et al. (1970); 
O'Neil et al. (1971); and Onishchenko et al. (1971), which was a phenomenologically based theory written 
down partly using physical arguments. Frieman et al. (1963) and Simon and Rosenbluth (1976), on the 
other hand, championed the Hopf scaling, and attempted amplitude expansions of the kind summarized 
above. Unfortunately, such expansions are plagued by critical-level singularities which unavoidably haunt 
the higher-order equations of the expansion, as noted in Sec. IV.A. To avoid this difficulty Frieman et al. 
(1963) and Simon and Rosenbluth (1976) removed the singularities by displacing the poles from the real 
axis in a regularization procedure that was meant to be equivalent to Landau's solution of the initial-value 
computation. But the procedure has questionable mathematical validity (strictly speaking, the asymptotics 
introduces multiple timescalcs to account for a relatively fast adjustment of the system to the slower growth 
of the unstable mode; using Landau damping ideas to remove the singularities in the expansion obscures 
the separation of timescales). 

Because neither approach was completely convincing, the effort to distinguish between them fell to nu- 
merical simulation and laboratory experiments, which suggested that trapping scaling dictates the saturation 
level (e.g. Denavit, 1981). Numerical computions of the bump-on-tail instability are presented in Figs. 4 
and 5. Figure 4 gives an illustration of the nonlinear dynamics occuring during the saturation of an unsta- 
ble mode. Figure 5 collects together results from a suite of computations that measure the scaling of the 
saturation level, and confirms trapping scaling. 
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FIG. 4 Vlasov simulation, a = 0.23 and k = 1 {L = 2tt). The first four pictures (left to right, top 
to bottom) show snapshots of the total distribution function, F[v) + f{x, v, t) as surfaces above the 
(x, f;)-plane at the times t = 140, 200, 400, and 1800. The penultimate picture shows the amplitude 
of ^{t) := {J^ dx 99^/-L)^/^; the dashed line shows the expected exponential growth of the unstable 
mode and the dotted line indicates the prediction of the single-wave model. The final panel shows 
the spatial average of + / at the beginning and end of the computation. 
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FIG. 5 Numerical data for the satura- 
tion measure (the first maximum in 
the amplitude of the potential, <&(t)), 
plotted against the distance to the sta- 
bility boundary, a — uq. The predic- 
tion of the single- wave model (based on 
trapping scaling) is shown by the solid 
line. 



Relatively recently, two articles have more convincingly demonstrated that the correct scaling is trapping 
(Crawford, 1995; del-Castillo-Negrete, 1998a, b). Crawford's approach follows the center-manifold route and 
attempts to derive an amplitude equation for the unstable mode. The approach begins with the Hopf 
scaling and runs into the technical complications associated with integrals that diverge at the wave-particle 
resonance as one approaches the point of neutral stability. Rather than treating the singularities with 
some ad hoc rule, Crawford rescaled the mode amplitude with sufficient powers of e in order to counter 
the divergent integrals. This procedure rescales the mode amplitude back to the trapping scaling and 
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avoids singularity, which formally establishes that trapping is the correct scaling. Unfortunately, Crawford's 
amplitude equation contains an infinite number of terms of comparable size and cannot be used to explore 
how the instability saturates. 

Del Castillo-Negrete opts for a different approach by transferring results from critical-layer theory of 
fluid mechanics (Churilov and Shukhman, 1987a) that exploit matched asymptotic expansions to heal the 
singularity. The scaling required in the matched asymptotics is equivalent to trapping scaling, and del 
Castillo-Negrete thereby derived the single-wave model systematically. This is also the tack we take below, 
but first we highlight the development of critical layer theory in fluid mechanics. 



C. Critical layer theory 

The 1960s also saw much development of weakly nonlinear theories in fluid mechanics following advances 
in describing the onset of fluid convection (e.g. Malkus and Veronis, 1958). The techniques were quickly 
adapted to viscous shear flow problem by Stuart (1960) and Watson (1960), but they could not be adapted 
for the corresponding inviscid problem because of the recurring critical-level singularities. Some progress 
was made by assuming that significant viscosity wiped out the singularity in a slender region surrounding 
the critical level (Churilov and Shukhman, 1987b; Huerre, 1987; Schade, 1964), an idea dating back to Lin 
(1945). These insights gave; some physical support to the idea that the poles in the asymptotic expansion 
could be healed by viscosity. However, viscous effects also inexorably diffuse away from the critical levels, 
widening the critical layer as time proceeds (Brown and Stewartson, 1978). Overall, this so-called viscous 
critical-layer theory glosses over the critical-layer dynamics and is not a suitable closer to the inviscid limit, 
the situation of relevance to this review. 

The appreciation that the dynamics could be richer led several authors to develop more complicated 
"critical-layer theories" . These theories all proceed by resolving the dynamics of the critical level in detail 
whilst simultaneously including various physical effects to relieve any divergences. For example, Benney 
and Bergeron (1969) resolved the critical-level singularity by adding nonlinear effects and thereby built 
steady nonlinear waves (some criticism of their detailed construction was given later by Haberman, 1973 
and Brown and Stewartson, 1978). Other theories were those of Stewartson (1978) and Warn and Warn 
(1978) for forced critical levels that incorporate unsteadiness and nonlinearity. In both examples, a key 
detail is that nonlinearity first becomes important where the solution is nearly singular, which is over the 
critical layer, the narrow region surrounding the critical level. 

To motivate the critical-layer theory we develop below, we again exploit numerical computations. Figure 
6 shows a simulation of instabilities growing on a shear flow in the form of a jet (the background flow is 
the so-called Bickley jet (Bicklcy, 1937), with U{y) = sech^y, which is popularly used to model the fluid 
wake behind an obstacle). Both this fluid mechanical example and the Vlasov simulation shown earlier 
show how the unstable mode grows initially, but then saturates as a result of pronounced nonlinear effects 
inside a narrow region surrounding the critical level. The principal effect in both examples is the twist up 
of the background distribution into vortices over the critical layer. Note how the perturbation vorticity 
distribution is clearly composed of two parts: a global, relatively smooth modal structure, and finely scaled 
cat's eye patterns emerging around the critical levels that coincide with the inflection points of U{y). 



V. THE SINGLE-WAVE MODEL 

The stage is now set for a reduction of our general model equations to the single-wave model. From 
Sec. Ill we learn that our linear theory has an eigenspectrum that is composed of only the continuum plus a 
single, distinguished, embedded neutral mode. Examples of such embedded modes include the special plasma 
oscillations lying on the stability boundary described detailed in Fig. 3. On adjusting a control parameter, 
and thereby perturbing the equilibrium profile slightly, that mode becomes weakly unstable. The single- 
wave model describes the ensuing nonlinear development, and an essential ingredient is the trapping scaling 
discussed in Sec. IV. 
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FIG. 6 Computed cat's eye patterns formed from the growth of an unstable mode in a simulation 
of the Bickley jet with k = 0.6. Computation is initializing with the third Fourier mode in x, 
and lower modes excited by round-off errors are explicitly removed. (Computational domain in y 
is much larger than indicated. Both x and y are periodic and a pseudo-spectral scheme with an 
explicit viscosity {u = 3.75 x 10^^) is used. ) Upper left is a plot of the background flow profile, 
and the vorticity, —dU/dy. Following (left to right, top to bottom) are five snapshots of 
the vorticity field at the times indicated. Comparison of the vorticity field snapshots with the 
—dU/dy plot indicates the shading scheme and contour levels. The bottom row of pictures shows 
the perturbation vorticity field, u, for three of the snapshots. 

A. Derivation 

At onset, a distinguished neutral mode is embedded within the continuum and remains regular because 
its critical level lines up with an extremal point of the equilibrium profile (where i?'(/*) = 0). At first 
sight, this might appear to allow us to proceed with a conventional weakly nonlinear theory as the critical- 
level singularity then disappears from the leading-order solution (the distinguished mode). However, the 
continuous spectrum manifests itself at higher order in the asymptotic expansion where the singularity 
reappears and becomes progressively more severe as one proceeds in the expansion (the core of the difficulty 
in early expansions for electrostatic instabilities exploiting the Hopf scaling of Sec. IV. B). The cure for 
the critical-level singularity is to abandon the conventional weakly nonlinear solution over a slender region 
surrounding the critical level. Within this critical layer, a different solution is needed that varies on a much 
finer spatial scale. The recipe is therefore a matched asymptotic expansion in which one finds an inner 
solution for the critical layer, and then matches to the usual weakly nonlinear solution that remains valid in 
the 'outer region' further afield. We sketch out this expansion for the general model; additional details can be 
found in the articles by Balmforth and Piccolo (2001); Churilov and Shukhman (1987a); del-Castillo-Negrete 
(1998a,b); and Goldstein and Leib (1988). 

To ease the discussion, we will take ^{I) to be monotonic, implying that there is only a single critical 



26 



level, and 17' > 0. Wc also assume that our periodic domain encompasses c;xactly one wavelength of the 
distinguished neutral mode and the system has been suitably scaled so that the wavenumber is A; = 1. In 
principle, the expansion requires multiple time scales in order to deal with the relatively fast propagation 
of the neutral mode as well as its much slower growth. However, if we first transform into the frame of 
the distinguished mode (9 — >■ 9 — <;t), where <: denotes the frequency of the frame shift, it then suffices to 
simply rescale and use the slow time, T = s~^t, where s is the small parameter that we use to organize 
the asymptotics. In other words, we begin by introducing the transformation, dt — *• —';d$+ sOt, into the 
governing equations (24) and (22), whereupon they become, respectively, 

epT + <;)ps - R'^e + [p, = , (75) 

$ = /Cp := / d/' / d(9' K{e - 9', I, I') p{6', I', t) ; 
Jv Jo 

these equations constitute the starting point for our expansion. 



1. Outer regular expansion 

Outside the critical layer we assume the asymptotic sequences, 

Ril) = Ro(I) + eRi(I) + . . . (76) 
p = e^p2{9, 1, T) + e''p3{9, 1,T) + ... (77) 
$ = /C/9 = s^MO, I, T) + e^'MO, I,T) + ..., (78) 

where i/jm = I^Pm for m = 2,3. The expansion in (78) opens at order which corresponds to the "trapping 
scaling" of O'Neil et al. (1971) discussed in Sec. IV. B, and the sequence for R{I) incorporates the distortion 
of the equilibrium profile needed to pass from the marginally stable state to an unstable one (such as by 
setting a = qq + e in the parameterized family of equilibria of the plasma problem of Sec III.C, as depicted 
in Fig. 3). 

Upon introducing the sequences into the governing equations and grouping together terms of the same 
order in e, we obtain to leading order e^, 

{n-^)deP2-Rode^2 = 0, (79) 
which can be integrated with respect to 9 to obtain 



^P2 



11 — c; 



P2 = 0. (80) 



The solution to (80) is the embedded neutral mode with Q(/*) = <j and R'o{I*) = 0: 

p^{9,I,T) = ^^, MO,I,T)=A{T)rj{I)e''+c.c., (81) 

where ri{I) is now the relevant (fc = 1) eigensolution to (39), and the complex amplitude, A{T), is not yet 
determined. Note that cr^^ri{I) is also the null vector of the operator adjoint to £. When evaluated at 
the critical level / = /« , the solution becomes 

P2{9,h,T) = ^<l>{9,T), ct>i9,T):=AiT)e'' + c.c., (82) 

where the subscript '*' notation refers to evaluation at / = /* and, for later convenience, we have normalized 
so that ??(/*) = 1. 

Proceeding to order we obtain 

(O - <;) depa - R'odeipa = Kd6^2 - dTP2 , (83) 
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which when integrating with respect to 6 and dividing through by — fl) yields 



fl-<; 

In general, although R'{I) vanishes at the critical level, neither R"{I) nor R'i{I) will vanish. Thus, the right- 
hand side of (84) is singular, which impedes us from multiplying by the adjoint e~'^?7(7) and integrating to 
eliminate pz to find the desired equation for A{T). To avoid this stumbling block, we resolve the critical- 
layer dynamics by using an inner solution. The scaling of the inner region is determined by noting that 
the asymptotic sequence, e^p2 + e^ps becomes disordered as I tends to J*. In fact, because p2 ~ 0(1) and 
P3 ~ 0{I — /»)^^ as / ^ the asymptotic solution breaks down within a slender layer of thickness e 
surrounding the critical level. This is where we must search for the inner, critical-layer solution. 



2. Inner critical layer solution 

For the critical layer, we define a stretched coordinate Y := {I — I^)/e and a new variable, ({9, Y, t), (not 
to be confused with the general C{q,p,t) of Sec. 11. A), such that 

p = s\{e, Y, T) + s^pT'^'^ie, h,T). (85) 

Because $ = ICp is an integral transform of p which smoothes that variable, there is no critical-layer structure 
in $ at leading order (cf. Sec. V.A.4), and therefore $ ^ e'^ip2(d, I*:T) = e'^4>{6,T) over this region. Given 
these variables, the relation 9/ = e^^dy, and Eq. (24), the equation for C, becomes, to leading order, 

St (C + ^ + Co - {R'u + Cy) h = 0, (86) 

where (A(^,T) is given by (82). 

Equation (86) is a fully nonlinear evolution equation for the critical-layer dynamics and offers no further 
simplification. The problem cannot be solved in closed form and must be approached numerically in general. 
Thus, the problem still boils down to the solution of a partial differential equation, and the single- wave model 
remains of infinite dimension. 



3. Matching 

Expressed in terms of Y, the inner limit of the outer solution can be written in the form, 



2 ^0*4' 2 ^1*4 2 ^0*4Te , 3 

^ 01 yoi Yfie 



+ 0{e') , (87) 



provided t/j^ = ICps does not diverge faster than (/ — /*) -'^ as / which is the case for the /C-operators 

of interest here. On the other hand, the inner solution has an outer limit given by 

^'C + e'^, (88) 



e 



with the second term being the outer contribution to the inner solution. The two leading-order solutions 
can therefore be successfully matched provided C for Y — > ±oo, which requires consideration the initial 
condition, ^(6*, y, 0). In particular, we take ({9,Y,0) = 0, which corresponds to kicking the system into 
action by introducing the discrete mode at a low amplitude, ^(0) = with |^o| ^ 1, without further 
rearranging the equilibrium distribution within or outside the critical layer. The large |y| limit of (86) then 
implies that 

Thus, C ^ 0{Y~^) for |y| 1, ensuring the match. 
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4. Solvability 

The final step is to apply a solvability condition on ps to determine A{T). To accomplish this task, first 
consider 

}Cp= [ d9' [ di' K{e-e',i,i')p{e',i',T). (90) 

Jo Jv 

We break the integral over V into three pieces: [II, — A], (J, — A, J* + A) and [/, + A, //{], where [Il,Ir] 
denotes the full integration over V and e <C A <C 1. The breakages in the integral therefore extend into the 
matching region where both the outer and inner solutions remain valid. For the integral in I over T> of (90), 
we then obtain 



£2 / di' K{0-e',i,i')p2{0\i',T) 

Jv 

+ eH dl' K{e-9',I,I')p3{e',I',T) 

J V 

.A/e 

+ e^K{0 - 0', dY C{0', Y, T) + O(e^) , 



where we now use the notation, 

n/.-A i-Ir 



I d7 /(/):= /* d//(/)+ r dlfil). 
J V Jil Ji,+a 



(91) 



Hence, 



^2 = / dl' K{e - 0', I, I') p2{0', I', T) (92) 
Jt> 



throughout V, which is therefore free of critical-layer structure as remarked earlier, and 

^2 

V'3 



= / d0 J- dl' K{0-0',I,I')psi9',I',T) (93) 

Jo J T> 

+ d0 K{6-0',I,h) dYC{e',Y,T). 

Jo J-A/e 

Now we multiply Eq. (84) by r]{I)e~^^ , integrate over [0, 27r], and / over the ranges — A] and 

[/* + A,7i{]. The left-hand side of (84) becomes, 

= -/ d0 dYe-'\{0',Y) + O{A), (95) 

Jo J-A/e 

on introducing tp^ from (93) and using the following properties of rj: C^rj = 0, rye"'^ = /C[iio77e~'^/(n — <,)] 
and r?* = 1. The right-hand side of (84), on the other hand, becomes 

Thence, 

= -— / d0 dre-'\(0,r,t) + o(A). 

Finally, we take the limits A — > and A/e oo. This turns the integrals over / in (97) into standard 
principal- value integrals, and the integral over Y into another principal- value integral, this time at its infinite 
limits. 
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B. Normal and degenerate forms 

Gathering together the results of Sec. V.A, we rewrite Eqs. (82), (86), and (97) compactly as 

(/.(^.T) =y^(T)e'^+^(^)*e-'^ (98) 
QT + [Q,e]=0, (99) 
imAT + i^A = -{Ce-''^) , (100) 



where 

Q(e,y,T) :=C+^0+^F' + <F, (101) 



D// D// 

_ _ 0* ^ , -"-0* 
( 



e{e,Y,T):="-fY^ + (P{e,T), (102) 

These equations can be simplified further by scaling, which we do in the next section to obtain the normal 
form. Before setting about this task, we briefly recall the origin of some of the quantities appearing in the 
system: the variable Q corresponds to the limit of the total distribution, C, = R(I) + p, the Hamiltonian 
system of Sec. II. A within the critical layer. The equilibrium distribution, R{I), contributes the coefficients, 
i?o« and (the leading-order curvature and the first-order slope), to the definition of Q. Featuring in 
(102) is the leading-order 'shear' at the critical level, f2'^ = h"{I^,), determined from the one-particle piece 
of the Hamiltonian. The integral quantities, w and p, depend on the shape of the eigenfunction of the 
embedded mode. 



1. Normal 

To place the system into a simpler form, we introduce a characteristic timescale r, a scale for action /x, 
and rescaling of the distribution, b. We then define the new variables, 

T 

t:=-, x = e + iJ.t, y:= n'^rY + n, (104) 

T 



A:=-n:TUe-"'\ (105) 
bC, Q = bQ- lix^K + nj, (106) 



where 



_ b_R^^ 



«^=-7Z7^> 7:=M«+^. (107) 



If Ci7 7^ 0, we may then make the selections, 



T V 

— , M := r- , (108) 

W ZD 



where w and v are given in (103); the "degenerate" case, with = will be dealt with separately in Sec. 
V.B.2. Substituting these forms into Eqs. (98), (99), and (100), and then dropping the on ^ leads to 

(/j(x, t) = yle'^ + A*e-'^ , (109) 

Qt + \QA = ^ (110) 
iAi = (Ce--), (111) 
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where 

Q{x,y,t):=C + K^-'^y'^+^y (112) 
£{x,y,t) := - ^{x,t) . (113) 

An alternative form of Eq. (110) is 

Ct +yCx + 'PxCy = -K.'Pt - I'Px ■ (114) 

The initial and boundary conditions are 

A{0)=Ao, C{x,y,0)=0, 

C ~ y-\Kip^t - IV) for \y\ > 1, (115) 

where Aq is the initial mode amplitude. Note that the relatively weak decay of C with y does not interfere 

with the construction of (Cc"'"^) because the integral over y is interpreted using the principle value at its 
limits (c/. (97)). Similarly, one can eliminate the variable Q entirely in favor of Q by observing that 

{Qe-'^ - kA) = {{Q - K^)e-'^) 

= ((C-f.^+7.)e-) 

= (Ce-'^) , (116) 

provided the potentially divergent integral involving ^Ky^e"'^ is dealt with by first performing the x-average. 

Equations (109) (113) are what wc call the normal form of the single-wave model. Crucially, the absence 
of dissipation in the problem has the unappealing consequence that no dimensional reduction occurs as 
part of the asymptotic analysis, and the single-wave model retains a partial differential form. Of course, no 
dimensional reduction can be expected here, since at the onset of instability, the distinguished neutral mode 
is arbitrarily close to an infinite number of singular continuum modes with comparable wavespeeds; the two 
equations composing the single-wave model describe the interaction between the distinguished mode and a 
superposition of locally resonant continuum modes within the critical layer. 

Note that one further simplification is possible: since r has not yet been specified, we may exploit this 
timescale to set 7 = ±1. The parameter 7 is built from the correction to the equilibrium profile, R\, and, 
as shown below, acts as the trigger for instability. After scaling, the two choices, 7 = ±1, refer to whether 
the system is displaced from marginality into the unstable (7 = —1) or stable (7 = +1) regimes. Aside 
from the initial amplitude ^0, this leaves a single parameter in the problem, namely k, which is related 
to i?Q^ , the curvature of the leading-order equilibrium profile at the critical level. This parameter encodes 
information regarding the position of the embedded mode along the stability boundary: at the minimum of 
the stability boundary where all wavenumbers are stable, i?o* = (see Sec. III.C) and so k = 0; elsewhere 
on the stability boundary, n docs not vanish (and the danger of additional unstable modes is removed by 
suitably choosing the length of the periodic domain). 

The normal form as given by (109), (111), and (114) was derived for various types of fluid shear flows 
by Balmforth (1999); Balmforth et al. (2001a); Balmforth and Piccolo (2001); Churilov and Shukhman 
(1987a); Goldstein and Leib (1988); and Shukhman (1989). The system is similar, though not identical, 
to the original single-wave model proposed for unstable electrostatic waves (O'Neil et al., 1971; Tennyson 
et al., 1994) and nonlinear Landau damping (Imamura et al., 1969). In the original single- wave model, the 
right-hand side of (110) did not appear (i.e., k = 7 = 0) and instability was introduced into the system via a 
suitable initial condition, C,{x, y, 0) ^ 0. O'Neil et al. (1971) used a beam-like distribution (a delta-function 
at a fixed velocity); Onishchenko et al. (1971) and Imamura et al. (1969) used a linear profile, equivalent 
to the K = case. del-Castillo-Negrete (1998a,b) derived the single-wave model above, but then elected 
to drop the parameters 7 and k and follow the route of O'Neil et al. (1971) to instability by selecting a 
suitable initial condition. However, when viewed as the weakly nonlinear description of the distinguished 
mode on the stability boundary, there is no freedom in the choice of initial condition; only the position on 
the stability boundary (and therefore k) is a free parameter. 
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2. Degenerate 

It may transpire that w = 0. When this happens, the scaling leading to the normal form of Sec. V.B.I 
is precluded because we cannot make the selections in (108). Upon examination of (103), it is evident that 
w = is unlikely unless there is a simplicity or symmetry to the problem, which turns out to be the case 
for some physical examples. A case in point is that of the XY model introduced in Sec. II. C. 3. For this 
case, ~ I and the interaction has a spatial kernel (— 7r~^ cos(d — 6') is independent of action /) so rj is 
constant (see Sec. III. A). Moreover, for equilibrium profiles like the Gaussian of Sec. III.D, R'q{I) is odd in 
/, the embedded mode has c = 0, and tu = follows immediately. We will treat this example in detail in 
Sec. VIII, while here we work out the general degenerate form. From Eq. (100) it is clear that the mode 
amplitude for the degenerate form will no longer be governed by a differential equation, but by an integral 
constraint. 

Instead of (108), we now make the selections 

»=^- "'-^ 

With these choices, 7 = and Eqs. (98), (99), and (100) now become 

^(a;,t) = Ae'=" + A*e-'^, (118) 

Qt + [Q,E]=Q (119) 

A=(Ce-'-), (120) 

where 

Q(a;,t/,t) :=C + K^- (121) 

£{x,y,t) := - ^{x,t) . (122) 

Finally, we are again free to choose r, and thereby scale out a parameter. The only parameter remaining 
explicitly is k, and so we choose 



1 

T = 



(123) 



so that K = — sgn(i?Q^/z/) = ±1. With this selection, as demonstrated in Sec. VI, we find that the system 
is unstable when k = —1 (and stable if k = 1). Note that, because the mode amplitude is no longer a 
dependent variable, but given by the integral constraint (120), the initial amplitude of the mode does not 
appear as an independent parameter, and the initial state of the system remains to be specified. One option 
is to initialize using the linear normal mode of Sec. VI. Alternatively, a more convenient method in practice 
is to kick the system into action by adding a transient forcing term to (120) that smoothly switches on the 
dynamics. 



C. Hamiltonian structure and conservation laws 



1. Normal Hamiltonian structure 

Both the normal and degenerate forms of the single-wave model of Sec. V.B inherit the structure of their 
Hamiltonian parent models, discussed in Sec. II. For the normal form of the model, the Hamiltonian is given 

by 

Hn= I dy I dx {£Q + K\Af^ . (124) 
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The two terms in (124) arc reminiscent of the particle kinetic energy and the electrostatic energy for the 
Vlasov-Poisson system (see Sec. II). This Hamiltonian together with the Poisson bracket of (7) gives the 
single-wave model. To see this we compute the derivatives, 

and ^ = -2n{Qe--.A), (125) 

where Q, A, and A* are independent variables, and insert them into (7) to obtain 

Qt = {Q,H^} = -[Q,£], 

At = {A, Hn] = -i (Qe-'^ - kA) (126) 

and similarly A^ ~ (At)* ■ This Hamiltonian structure is the continuum version of that derived in Tennyson 
et al. (1994) (see also Mynick and Kaufman, 1978; Smith and Pereira, 1978). 

Showing that (124) exists requires some care, of the type commented upon in Sec. V.B.I, since the second 
term of Hn divergent. To show that Hamiltonian (energy) is in fact finite we return to (114) and 
analyze more carefully the large \y\ behavior of (, obtaining 



C ~ ^{^Vxt - IV) + -^{lA-^? - KiPtt - IVxt) + - ■ (127) 

y 



The a;-independent, 0{y^^) term in (127) follows upon averaging (114) over x, and is sufScient to cancel 
the divergent second piece of Hn, leaving a convergent double integral. 

Given the Hamiltonian structure, it is natural to look for invariants. Since both the single wave and 
the medium can carry momentum we expect a total momentum of the following form to be conserved (cf. 
Tennyson et al. (1994)): 

P = 2Tr\A\^ - [ dy [ dxyC, (128) 
Jr Jo 

where the first term of (128) will shortly be identified with the wave action, while the second has the familiar 
form of the (relative) x-component of the momentum in vortex dynamics (e.g. Flierl and Morrison, 2011). 
It can be shown that P commutes with H^, i.e. {P, H^} = 0, or directly that this is equivalent to 

Pt = {{yO-\A\')^=0. 

An immediate consequence of the Poisson bracket of (7) is the existence of the following infinite family 
of Casimir invariants (see e.g. Morrison, 1998): 

C = [ dy I dxC{Q) , (129) 



where C is an arbitrary function of Q. An important Casimir is the 'enstrophy' 



2tt 



C2 := dy dxQ\ (130) 



With some manipulation this can be shown to be equivalent to 



C2= dy dxC"^ -2KHn + 2j{2Tr\Af - P) . (131) 

Jo 



Since and P are constants of motion, this also implies the further invariant, 

dy / dxC^ +4Trj\Af . (132) 
Jo 

From (132) we obtain a stability condition, viz., for 7 — +1 in the vicinity of the equilibrium state A = C, = 0, 
the invariant / provides a norm for bounding solutions to an infinitesimal neighborhood of this equilibrium. 
Also, with the initial conditions ({x,y,0) = and ^(0) = Aq we have 

7|A|2 + ^(a=7l^or (133) 
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Thus, for 7 = +1, wc have the finite amplitude bounds [^p < \Aq\'^ and (C^) < 2|j4op. This amounts to 
energy-Casmir stability for this example (see e.g. Morrison, 1998). 

To show that J = 2tt |Ap is equivalent to the wave action, we introduce the coordinate change. 



2tt 



and 



(134) 



The bracket of (7) becomes 



{F,G} 



dFdG _ dGdF\ 

dy dxQ 



SF SG 
SQ'SQ 



(135) 



where now \1/ and are a genuine canonically conjugate pair of wave variables for the wave. In terms of 
these variables the potential is given by 



ip = 



cos{x + 4*) 



and Eq. (Ill) is replaced by the following two real-valued equations: 



* = 



J 



dHn 
dj ~ 

_dHn 
5* 



■ Q cos{x + ^) — K 



= -2 Q sin(a; + *) ) . 



(136) 



2. Degenerate Hamiltonian structure 

It is evident from Sec. V.C.I that divergences associated with system size must be treated with some 
delicacy; e.g., observing that (</?) vanishes requires integration over x, giving zero because of periodicity, 
before integration over which woiild give a divergent result, as discussed in obtaining Eq. (116). Such 
system size divergences are common for two-dimensional theories such as vortex dynamics, quasigeostropy, 
etc., which possess infinite kinetic energy for unbounded domains. For this reason, relative energies and 
momenta arc often used (sec, e.g., Flicrl and Morrison, 2011 for an example of numerically tracking such 
quantities), but for our presentation of the Hamiltonian form for the degenerate case, we find it convenient 
to introduce a system size so that (1) = T>. The need for this artifice is compounded by the fact that 
for the degenerate form, the theory is written entirely in terms of the variable Q - the amplitude being 
determined by Eq. (120), which is a constraint equation. Tools for handling such constraints are now well- 
developed (Chandre et al., 2013, 2012; Tassi et al., 2009), but we will present the Hamiltonian structure for 
this degenerate case in a manner akin to that for the Vlasov-Poisson system (Morrison, 1980), which has a 
bracket of the form of (12). 

Using (118) and (120) it is evident that the constraint amounts to if being functionally dependent on 

i.e., 

ip{x) = 2(Ccos(a; - x')) =: ip[C]{x) (137) 
where the underline will be used to denote if as an operator. Upon introducing the linear operator 

Ci^ :=I + Kip, 

where I is the identity operator, the equation of motion (119) becomes 

C4Ct] = -yCx-l[C]xCy (138) 

and thus it is necessary to assume £„[Ct]^^ exists, which we write as 



r — l l -7- I 2 2 

= = X — KW + K W — 



(139) 
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Various properties of the operators Lp and will be needed, and so we present them now. Observe, 
both operators are formally self-adjoint, e.g, 

{g^h\) = {h^[g]) 
for functions g,h{x,y), and we have the following: 

£ = 2/V2 - ^[C] 



(140) 
(141) 
(142) 



Also, by periodicity, (p[l] = (fi[g{y)] = 0, where 5 is a function of y alone. A short calculation gives 

^ [^[C]] = v^[c] ^"[C] = V[C] • (143) 

Thus, upon making use of (139) we obtain 



and a functional chain rule calculation (e.g. Morrison, 2013) yields 



5F 



5Q 



6£ 



(144) 



(145) 



which we will find convenient. 

Now, the Hamiltonian of (124) for the normal form suggests a Hamiltonian for the degenerate case with 
terms as follows: 



(1 + kV) 
2tt 



dx y C 



dy dx dy' dx' 
V Jo Jv Jo 



X C(a;, y)C{x, y') cos{x - x ) , 
which has a form reminiscent of the Hamiltonian for the Vlasov-Poisson system. Evidently, 

5Hd y" iX^nV) 



(146) 



'-^ [ dy'f dx'C{x',y') cos{x-x') 

TT J-D Jo 



(1 + kVMC] 



Making use of (147), (145), and (147) gives 



SJh^y^ 
SQ 2 



whence, we obtain 



Qt = {Q,Hd}^-[Q,£] 



(147) 



(148) 



(149) 



the equation of motion of (119). This Hamiltonian structure is similar to that given in Bachelard et al. 
(2009) for the XY model. 
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D. Dissipative effects witfiin the critical layer 

Depending on the physics under consideration, the Hamiltonian forms described in Sec. V.C can be 

amended by adding terms that describe dissipative effects. In particular, dissipation terms can be included 
on the right-hand side of the critical-layer equation (114), corresponding to various kinetic theory transport 
processes: 

Ct + yCx + ^xCy + Kipt + -yifix = (ACy )!/ • (150) 

The dissipative term corresponds to viscous dissipation in the fluids context, or a Fokker-Planck operator 
in the kinetic theory context. For the latter, A may depend on the coordinate y. 

The inclusion of such dissipative terms within the ordering scheme of Sec. V.A requires justification. In 
each case, there needs to be a distinguished scaling that ensures the term enters at the same order as the 
other physical effects inside the critical layer, whilst simultaneously not affecting the outer solution. We do 
not pursue this here in generality, but for the fluid problem the kinematic viscosity v of the parent model, 
which provides a term viuixx + ^^yy) on the right-hand side of Eq. (28), must be scaled as u = e^u^, and A 
is then a suitably scaled version of (cf. Churilov and Shukhman, 1987a; Goldstein and Hultgren, 1988). 
Note that if A ^ 1, one can perform some additional reductions of the model to furnish simpler amplitude 
equations (see Churilov and Shukhman, 1987a), but here we focus on the less viscous cases with A order 
unity or smaller. 



VI. LINEAR SINGLE-WAVE DYNAMICS 

We now explore the linear theory of the single-wave model, focussing mostly on the normal form of Sec. 
V.B.I. On discarding the nonlinear term, fxCy^ the system can be written as 

Ct+yCx + Kipt + JiPx = (ACy)y , (151) 

lAt = (C e"'^) , (fi = Ae'^ + c.c. , (152) 

where we have included the dissipative terms of Sec. V.D. These equations can be solved analytically for 
both normal modes and the initial-value problem. As a result of the latter, we arrive at one of the simplest 
examples of a system exhibiting Landau damping. 



A. Normal modes 

We consider first the ideal case where the dissipative terms of Sec. V.D are set to zero (i.e. A = 0), and 
then we investigate what happens when we include them. 



1. Ideal 

Assuming that nondissipative normal modes have the dependencies, cx e^i^-'^^-*) and A cx e^"^^' (recall 
our single-wave has wavenumber k = 1), where ci is the (complex) wavespeed, Eqs. (151) and (152) reduce 
to 

= '^^i_Zj ^e^^ + C.C. (153) 

y-ci 

and 

c\A = (kci — 'y)A [ dy — - — = i7r(K;ci — 7)^ sgn(cij) . (154) 

Jr y-ci 

The second equality of (154) is valid for all Cn ^ 0; it is easily obtained by closing the contour with a 
semicircle in the upper or lower half plane, depending on the sign of Cu, and using Cauchy's theorem. 
Solving (154) for ci gives 

•^1 = 1 ,^^2 2 [^'^ ~ i sgn(cH)] , (155) 



36 



Thus, cij = — 7r7 sgn(ci,;)/(l + 7r^K^), which is inconsistent unless 7 < 0. Since we have normahzed 7 = ±1, 
this indicates that the system is unstable when 7 = —1. However, if 7 = +1, there is no consistent solution 
for cii, and so there are no normal modes and the system is stable. 

The factor of sgn(cij) in (154) arises because the integral has a branch cut along the real axis of the 
complex ci-plane, which as noted in Sec. Ill reflects the continuous spectrum lying along the contour of 
integration. If we start from the consistent unstable solution with cu > (7 = —1) we can follow Landau 
and analytically continue ci into the lower have plane by deforming the contour of integration so as to 
capture the pole. This gives 



dz = in 

Z — Cl 



where the symbol on the integral is selected to show the contour is deformed so as to enclose the descending 
pole. This procedure of Landau amounts to analytically continuing through the branch cut and onto a 
different Riemann sheet with cn < 0. Thus, the analytic continuation of the dispersion relation in (155), 
which would appear in the long time limit of the Laplace transform solution of the initial- value problem (as 
in Landau's classical analysis), has no sgn(cij). The damped solutions one obtains in this case for 7 > 
correspond to the quasi-modes discussed in Sec. III.C. 

Following suite for the degenerate normal form of Sec. V.B.2, we find that the normal-mode solutions 
satisfy 

y-ci 

and 

f 1 

A = KCiA / dy = iTTKCiA sgn(cii) . (157) 

Jr y - Cl 



That is, 



c = \Ci, \ci\ = — (158) 

■KK 



Hence, there are normal modes and the system is unstable if k = — 1. 



2. Dissipative 

Assuming that A is constant, the dissipative modes satisfy the following equation: 

- i(y - ci)C = 1(7 - «ci)(p . (159) 

Equation (159) can be solved by using Airy functions or Fourier transforms (generalizing the case with only 
the viscosity term as in Benney and Bergeron, 1969 and Balmforth, 1998): 

C = i(KCi - 7)<^ J{y, ci; Ao, Ai, A2), (160) 

where 



1*00 

J:= dge-^«'/3-K-i«(y-ci) 
Jo 



Upon inserting (160) into (152), and then interchanging the order of integration, we obtain S{q) for the 
integral over y. Using this to evaluate the integral over q gives the eigenvalue, 

Cl = ^ ""^^ . {ttk - i) ■ (161) 

Equation (161) is identical to the inviscid result, save that sgn(cii) no longer appears. Thus, if 7 = —1, 
the system is unstable, and the inviscid and viscous growth rates coincide. But when 7 = -|-1, there is 
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FIG. 7 Eigenfunctions for k = 0. Panel (a) shows the unstable mode for A = 0.1 and 0, and panel 
(b) the decaying modes for A = 1, 0.1 and 0.01. Also shown in (b) is the asymptotic result of 
Eq. (162). 



a decaying viscous normal mode with no inviscid counterpart. However, it does correspond to the fake 
eigenvalue (quasi-mode) obtained by analytical continuation of the inviscid dispersion relation. 

The eigenfunctions associated with the unstable and decaying eigenvalues are illustrated in Fig. 7. The 
unstable viscous mode has a regular shape and smoothly limits to the inviscid eigenfunction. However, 
the decaying viscous modes oscillate spatially with a wavelength that becomes ever shorter as the viscous 
parameter A is decreased. For A ^ 1, the approximate shape of the decaying eigenfunctions can be found 
using the method of stationary phase on the integral in (160): 



^/^T{KCl - j)(p 



exp 



5i7r 



Al/4(y_ci)l/4 

Clearly, this viscous mode has no regular inviscid limit. 



(162) 



B. The initial-value problem and Landau damping 

The solution of the initial-value problem for C can be written in the form, 

rt 

( = -e^^ / [kA,{s) + i7A(s)]e-'^(*-")-^(*-'*^'/3 + ^ ^ (163) 
Jo 

Again, the integral over y generates a delta-function, allowing the mode amplitude to be found, 

A = Aoe", r = .. ^. , (164) 

(iTTK — 1) 

where F is the viscous normal-mode eigenvalue. The complete solution is then given by 

Aoer*+'^ n 



C = ° , / e-r«-'^«-^^ /Ms + c.c.. 

ITTK -1 Jo 



(165) 
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-1 1 

X 

FIG. 8 Snapshots of the solution to the hnear version of the single- wave model at the times, t = 0.1, 
0.4, 0.75, and 1.5, showing the tilting of a initially periodic disturbance in a stable flow (7 = 1 
and K = 0). The right-hand picture shows the increasingly crenellated (phase mixed) distribution 
along X = 0. 



For A — > 0, Eq. (165) reduces to 

i^^oe'"(e"-e-i^*) , 

C = : \- c.c, (166) 

7r7 — TTKy — ly 

Alternatively, for i S> 1, we may replace the upper limit of the integral in (165) by infinity to leading order, 
recovering the viscous eigenfunction in (160). Thus, if the system is unstable, the solution of the initial- 
value problem converges to the inviscid or viscous normal mode. When the system is stable, on the other 
hand, the initial-value solution converges to an exponentially decaying disturbance which is equivalent to 
the decaying viscous mode. In the inviscid problem, the decaying disturbance corresponds to a Landau pole 
or quasi-mode. It is clear from (166) that the quasi-mode has a form which is not separable in space and 
time (the combination i{x — yt) appears), and consequently it cannot be a true normal mode. Moreover, 
the quasi-mode is evidently the inviscid limit of the decaying mode, as found in some related problems 
(Balmforth, 1998; Stewartson, 1981). 

Mathematically, Landau damping results from phase mixing in integral superpositions of modes of the 
continuous spectrum. The solution in (164) and (166) is, in fact, one such superposition. In (166), we observe 
the non-decaying, non-separable factor, expi(a; — yt), reflective of the basic advection and the continuous 
spectrum (see Sec. III. A). This non-separable structure of the superposition is also directly responsible for 
Landau damping: the long-time limit of integrals of the form, 

/ AyF{y)e^^-^-y'\ 
Jr 

vanishes according to the Riemann-Lebesgue lemma (Stein and Weiss, 1971), provided |i^(y)| is integrable. 
Moreover, in the case where F{y) is analytic in a strip containing the real y axis (as is the case for (166)), 
the decay is exponential, and this the mathematical essence of Landau damping. 

Physically, the non-separable superposition of continuum modes describes how a perturbation in the 
plasma or fluid tilts over as a result of the action of particles orbiting at different speeds or from the 
shearing flow of fluid elements. A tilting, initially periodic, disturbance then evolves into an increasingly 
crenellated (phase mixed) structure as shown in Fig. 8. There is no decay in the tilting process: the plasma 
distribution function or fluid vorticity is simply re-arranged. However, integral averages of these quantities, 
such as the electric field or streamfunction, do decay as a result of increasing cancellations. This does not 
imply that energy dissipation takes places; merely that energy is transfered between the various components 
of the energy as shown in Balmforth and Morrison (2001) and Morrison and Pfirsch (1992). 

Note that, although exponential decay dominates the decay of the streamfunction in the single- wave 
model, shear tilting normally leads to algrabraic decay in fluid shear flows (e.g. Brown and Stewartson, 
1980). A full solution of the linear initial- value problem for shear flow does indeed uncover both algebraic 
and exponential behavior (Briggs et al., 1970). The former disappears in the single- wave expansion because 
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it lies at higher order in the asymptotic scheme. Thus, smaU disturbances to stable flows would decay 
exponentially for a period, as the quasi-modes appear and subside, but ultimately a slower, more protracted 
algebraic decay sets in. This behavior is seen both in solutions of the linear initial-value problem for 
difl'erentially rotating vortices, and in electron plasma experiments designed to imitate them (Schecter et al., 
2000). 



VII. NORMAL FORM PATTERN FORMATION 

We now consider the nonlinear patterns described by the single-wave model; we first consider the unstable 
case (7 = —1), and then discuss the stable case (7 = 1). 



A. Pattern formation via unstable modes 

A typical example of the nonlinear evolution of an unstable mode is shown in Fig. 9, computed numerically 
(for a description of the scheme, see Balmforth et at, 2001a and Balmforth and Piccolo, 2001). Initially, A{t) 
grows exponentially as predicted by linear theory. However, when the mode amplitude reaches order one 
values, that growth is interrupted and the linear instability saturates. The saturation occurs when the mode 
twists up the total critical-layer distribution Q of (112) into a localized pattern. Beyond saturation, the mode 
continues to wind up the distribution, creating a well-defined cats-eye structure in Q that corresponds to 
vortices in the fluid and electron holes in the plasma. The flnal x-averaged Q in Fig. 9 shows a suppression of 
the original equilibrium gradient at the center of the critical region. Thus, saturation occurs when the mode 
removes the unstable mean gradients surrounding its critical level. Importantly, the flattened equilibrium 
distribution is contained well within the critical layer, and therefore occupies a much narrower region that 
the area spanned by the imstable gradients of the equilibrium (e.g. the bump in the plasma problem). A 
complete suppression of those gradients requires much more extensive re-arrangements, which could only be 
achieved if the mode amplitude saturated at a level dictated by the Hopf scaling. 

Overall, the dynamics culminating in the formation of the localized structures matches that seen with 
full simulations. The electron hole formation displayed in Fig. 4 is ubiquitous in Vlasov-Poisson simulations 
of unstable equilibria (see also, e.g., Cheng and Knorr, 1976; Heath et al., 2012 and references therein). 
Similarly, the cat's eye structures seen in Fig. 6 are a common pattern that forms in simulations of un- 
stable shear flow. A benefit of the single-wave model is the attainment of finer resolution with equivalent 
computational cost. 

Saturation also does not produce a steady mode amplitude. In fact, aperiodic "trapping" oscillations 
ensue (as they are often called by plasma physicists, since they reflect oscillations in the population of 
particles that are resonantly trapped in the potential created by the wave). The initial trapping oscillations 
are connected to the first overturns of the level sets of Q. These oscillations subside as the localized 
structure (the cat's eye) becomes well developed. The persistent trapping oscillations that occur much 
later appear to correspond to the sporadic formation of secondary structures within the localized structure 
that are sometimes called "macro-particles" (del-Castillo-Negrete and Firpo, 2002; Tennyson et at, 1994). 
Unfortunately, numerical computations with finite resolution are unable to determine the ultimate fate of 
the system and its trapping oscillations. Phase mixing ideas (cf. Goldstein and Leib, 1988; Killworth and 
Mclntyre, 1985; O'Neil, 1965) suggest that the secondary structures eventually become sheared out and 
disappear to leave a smooth, steady cat's eye. However, numerical simulations do not show much sign of 
the expected algebraic decay and numerical diffusion obscures the true dynamics (Balmforth and Piccolo, 
2001). 

In fact, explicit dissipation like the viscous term in (150) has an especially destructive effect on the cat's 
pattern: when one adds viscosity, the rearrangements of the distribution within the critical region begin to 
diffuse outwards. Consequently, the cat's eyes spread and the critical layer thickens continually with time. 
This allows the mode to access fresh, unstable equilibrium gradients, prompting a secular growth in A{t); 
asymptotic estimates and numerical computations indicate that A ~ t^/^ over long times (Balmforth and 
Piccolo, 2001; Churilov and Shukhman, 1987a; Goldstein and Hultgren, 1988). Secular growth continues 
until the system passes out of the the single-wave regime. In this situation, the single-wave model only 
captures an initial transient and the trapping scaling is eventually broken. 
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FIG. 9 Computation of an unstable mode in the single-wave model for 7 = — 1 and k = —0.4. (a) 
shows a time evolution of Ar (dashed), Ai (dotted) and ±|^| (solid), and (b) shows the initial and 
final x-averages of the total critical-layer distribution, Q = —Ky^/2 + jy + Kip + (. The stars in (a) 
denote the times at which snapshots of q are shown in (c) as densities on the (x, y)-plane, 

B. Pattern formation in stable systems; nonlinear quasi-modes 

In dissipative systems, the occurrence of instability and the resulting dynamics is usually the most 
interesting problem; stable systems relax back to equilibrium exponentially quickly and the only complication 
that can occur is if there are multiple possible end-states that "frustrate" the evolution of the system. 
Nondissipative systems like the pendulum, on the other hand, have an infinite number of possible periodic 
orbits adjacent to the stable equilibrium position, a general property of finite-dimensional Hamiltonian 
systems (see, e.g., Moser, 1976), and the Darwinian choice amongst them is determined by the energy of 
the initial condition. Stable ideal plasmas and inviscid shears also support an infinite number of adjacent 
wave-like nonlinear states. For Vlasov plasmas, the explicit construction of such states (the so-called BGK 
modes) dates back to Bernstein et al. (1958); in shear flows, Kelvin-Stuart cat's eye patterns (Drazin 
and Reid, 1981) are other examples. Nevertheless, despite this wealth of equilibria and the absence of 
dissipation, linear theory predicts that stable ideal plasma and inviscid shear flow relax through Landau 
damping. However, using the single-wave model, we now illustrate how sufficiently nonlinear perturbations 
do not continually decay. 

Numerical solutions of the nonlinear initial- value problem with 7 = +1 and k = and differing initial 
amplitudes are shown in Fig. 10. When the initial amplitude is sufficiently low, the numerical solutions do in- 
deed display increasing crenellation of the distribution and the mode amplitude decays as predicted by linear 
theory. However, for larger initial amplitudes, the damping becomes interrupted, and the mode amplitude 
returns to larger values where it fluctuates unsteadily in a sequence of aperiodic trapping oscillations. 

The interruption of Landau damping corresponds to the creation of another cat's eye pattern. Thus, 
sufficiently strong nonlinear perturbations realize Kelvin's objections (sec Sec. III.B) to inviscid shear flow 
dynamics. In the plasma context, the interruption of Landau damping has been observed in many compu- 
tations with the Vlasov equation and in experiments (e.g., Brunetti et al, 2000; Danielson et al., 2004; Feix 
et al., 1994; Franklin et al., 1972; Lancelotti and Doming, 2003; Manfrcdi, 1997; Sugihara and Kamimura, 
1972), and the nonlinear structures are assumed to converge to steady BGK modes. However, as for the 
unstable problem, the trapping oscillations persist throughout the computations in Fig. 10, and the ultimate 
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FIG. 10 Computation of a stable mode in the single-wave model for 7 = +1 and k = 0. Panel 
(a) shows scaled mode amplitudes, A{t)/AQ, for perturbations with different initial amplitude, Aq. 
The dots show the trend of Landau damping. Panel (b) shows the evolving cat's eye pattern for 
^4(0) = 7 at the times indicated by the circles in (a). In the particular case k = 0, the single- wave 
model possesses solutions with a real mode amplitude, A{t), and the reflection symmetry, y — )• —y, 
X — )• —X, and C. — )• — C- 



fate of the system is again unclear. 

The critical threshold in Aq below which Landau damping continues unabated occurs for initial am- 
plitudes of about unity. However, this threshold does not provide a clean criterion for the formation of a 
cat's eye because for slightly larger initial amplitudes (e.g. Aq = 2 and 3 in Fig. 10), the mode amplitude 
repeatedly switches sign. This implies that the mode twists up the distribution in one direction for a period, 
but then the sense of overturning switches and the mode unwinds the twisted background; it is not clear 
whether a coherent cat's eye ever forms. Only for somewhat larger initial amplitudes {Aq ~ 5 and 7 in 
Fig. 10) does the mode amplitude remain of one sign and the critical-layer distribution winds up continually 
into a cat's eye. 

With a simple reinterpretation of variables, the single-wave model can be applied to problems in polar 
coordinates (cf. Sec. II.C.l). Thus, thresholds observed in the single-wave model apply to the formation of 
nonlinear structures in desymmetrized vortices, both in numerical computations (Moller and Montgomery, 
1999; Rossi et at, 1997) and in non-neutral plasma experiments (DriscoU and Fine, 1990). The interest 
in these thresholds was sparked off by the idea that Landau damping (or, more generally, "continuum" 
damping, which includes algebraically decaying perturbations) could act to cause disturbed vortices to return 
to axisymmetry, rationalizing the preponderance of almost circular vortices in two-dimensional turbulence 
and atmospheric weather systems. However, if the vortex is too distorted, then it does not return to 
axisymmetry, a cat's eye develops and a multipolar vortex forms. Two numerical simulations of the two- 
dimensional Euler equations, showing the evolution of an elliptically distorted, Gaussian vortex are shown 
in Fig. 11. The more weakly distorted vortex displays sub-threshold axisymmetrization, whereas the larger 
distortion creates a multipole. By modeling the disturbance as the response of a weakly nonlinear quasi- 
mode and then deriving an appropriate version of the single- wave model to describe the dynamics, Balmforth 
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FIG. 11 Simulation results of elliptical perturbations to a stable (Gaussian) vortex. Shown are 
contours of constant vorticity. In panel (a) the perturbation is small and the vortex axisymmetrizes, 
while in panel (b) the perturbation is above a threshold for the formation of a multipole state. 

et al. (2001a) predicted that the critical amplitude threshold translates to an initial aspect ratio of about 
1.1 (a ten percent departure from circularity); the two simulations in Fig. 11 do indeed straddle this value. 
Turner and Gilbert (2007, 2008) present a much more thorough analysis of the threshold for desymmetrized 
vortices. 



VIII. XY MODEL; PATTERNS OF DEGENERATE FORM 

As noted in Sec. V.B.2 the XY model (also called the Hamiltonian mean field model by Campa et al. 
(2009) and others) provides an example of the degenerate form of the single-wave model; in this section we 
explore the pattern formation described by both the XY model and its degenerate single-wave form near 
onset. For the task, we consider the Gaussian family of equilibrium profiles, 

i? = i?o + £i?i = ^-^7^exp(-i/2) . (j^gy^) 
V 27r 

For this profile, — 0, R'^^ = 7 = 0. Moreover, fi'^ = 1 and (because this example has a spatial kernel) 
7] = 1, indicating 




(168) 
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The degenerate form of the single- wave model for this model is therefore 

Ct + yCx + VxCy = -i^Vt , A = (e-'^C) , (169) 

where k = — sgn(i//i?Q^). 

In Sec. VI, we established that the system is unstable when k = —1, and the growth rate is tt"^, or 
— r2'^i^/(7ri?Q^) — > (a — l)-\/2/7r, in terms of the original variables. When k = +1, the system is linearly 
stable and small perturbations decay through Landau damping; again, the exponent is (a — \)\f2fK for the 
original variables. Note that, except for the two choices, k = ±1, determining the stability of the system, 
the degenerate form of the single-wave model of (169) has no free parameters. 



A. Unstable modes and degenerate phase-space patterns 

Consider now patterns that occur in the XY model with the Gaussian equilibrium profile by increasing 
the amplitude parameter a through unity to trigger instability. We first solve the XY model numerically 
by using the operator splitting scheme of Cheng and Knorr (1976). We begin the simulation by kicking 
the equilibrium state by adding the forcing, 10"''fe^^°* cos6', to the potential <& in (33). This procedure 
prepares an initial condition for the dynamics after the forcing term decays to zero. Such prepared initial 
conditions, were called dynamically accessible in Morrison and Pfirsch (1989, 1990, 1992), where they were 
advocated because they preserve phase space constraints. 

Figure 12 shows results for such a kick-started simulation with a = 1.2. Specifically, Fig. 12(a) displays 
the temporal behavior of the amplitude, A{t), of the potential, defined in 

$ = A{t)e'^ + ex., (170) 

along with the behavior of linear growth (dashed line), for the Gaussian equilibrium depicted in Fig. 12(b). 
Figure 12(a) shows the solution initially following the linear growth, followed by oscillations in A{t) that 
are symptomatic of the twist-up of p{9, 1, t) into a cat's eye pattern around the peak of the distribution. In 
Fig. 12(c) we see that the corresponding average in 9 of the total distribution R{I) + p{9,I,t) is thereby 
flattened, suggesting how the instability is able to saturate itself nonlinearly. Not surprisingly, beyond 
saturation, trapping oscillations continue in A{t). Figure 12(d) depicts snapshots of the total distribution 
at the times indicated by the stars in Fig. 12(a), and shows how cat's eyes develop during saturation. 

Figure 13(a) shows scaling data obtained from a suite of computations of the type depicted in Fig. 12; 
the first maximum of A{t) is plotted (cf. Fig. 12(a)) against the control parameter a. This figure verifies 
the trapping scaling of the saturation level. 

The remaining two panels of Figure 13 show the numerical solution of the degenerate form of the single- 
wave model in (118)-(122). The system is again kicked into action by suitably forcing the potential (p. 
The unstable mode grows and saturates in a similar manner to the dynamics seen for the full XY model 
in Fig. 12. The mode amplitude, A{t), peaks first for Am ~ 0.52, which leads to a prediction for the first 
maximum of A{t) in terms of the original variables: Amax ~ 2TrAm{a — 1)^. This prediction is included 
in Fig. 13(a); the scaling data for the amplitude maxima from the full XY model computations converge 
satisfyingly to the prediction. Note that trapping scaling means that the rearrangment of p{9, 1, t) is unable 
to suppress the distribution to the point where the ^-average shown in Fig. 12(c) corresponds to a linearly 
stable equilibrium. 



B. Patterns in stable systems; nonlinear Landau damping in the XY model 

The parallels between the XY model and the Vlasov equation, together with the fact that both can be 
reduced to the single-wave model in the vicinity of onset, suggests that the XY dynamics should also display 
pattern formation even when a < 1, provided the system is forced with sufficient strength. That is, the XY 
model should display a nonlinear arrest of Landau damping akin to that discussed in Sec. VII. B. To confirm 
this expectation we once more solve the initial-value problem for the full XY model, again initiating the 
computations by applying a suitable kick to the potential. 
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FIG. 12 Kicked (dynamically accessible) solution of the XY model of (169) for the Gaussian 
equilibrium distribution of (50) with a = 1.2. Panel (a) shows the mode amplitude, A{t), along 
with the expected linear growth (dashed line). Panel (b) shows the equilibrium distribution as a 
surface over the {9, 1)-plane. Panel (c) shows the initial and final ^-averages of the total density. 
The stars in (a) indicate the times of the snapshots of the total density R{I) + p{9, 1, t) shown in 
the lower row as surfaces over the (0, /)-plane, focussing on the central maximum. 




FIG. 13 Numerical results for the XY model with Gaussian equilibrium and kick start (cf. Fig. 12). 
Panel (a) shows scaling data for the XY model, with the height of the first maximum in A{t), the 
amplitude of the potential, plotted against the control parameter a — 1. Also shown is the trapping 
scaling prediction of the single- wave model (SWM, dashed line). Panels (b) and (c) show the 
numerical solution of the degenerate single-wave model: panel (b) shows the mode amplitude, 
A{t), vs. time, t, along with the expected linear growth (dashed line); the stars of this panel mark 
times for the snapshots of the solution, Q = + ^ — (/?, shown as densities on the {6, /)-plane in 
panel (c). 



Figure 14 shows a computation in which a stable Gaussian equilibrium with a = 0.9 is forced by adding 
to $(f?,i) the external potential, 

*e = 2Aie-^°*'cos6l. (171) 
In Fig. 14(a) a plot of the logarithm oi A ^rs t is shown for a drive amplitude Aq — 0.5. This plot shows 
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FIG. 14 Solutions to the nonlinear XY model for a = 0.9, initiated by adding the external forcing 
(171) to the potential ^{9,t). In panel (a) = 0.5, and the stars indicates the times at which 
snapshots of p{9,I,t) are shown below in panel (c) as densities on the (^,/)-plane. In panel (b), 
the smaller values of = 0.1 and 0.02 are used, and longer intervals of Landau damping are 
observed. The dotted lines in (a) and (b) indicate the solution to the linear problem. 



how the initial decay becomes interrupted for a sufficiently large initial perturbation, as in Vlasov-Poisson 
simulations (e.g., Cheng and Knorr, 1976; Heath et at, 2012). Fig. 14(c) displays plots of the total density 
at the times indicated by the stars in Fig. 14(a); the interruption of Landau damping is caused by trapping 
and the emergence of the nonlinear quasi-modes visible in these images. These cat's eye structures are 
symmetrically displaced about the maximum of the original equilibrium profile at / = 0. Reducing the 
amplitude of the initial kick as in Fig. 14(b) prolongs the transient in which the system decays through 
Landau damping, but we have not attempted to compute any threshold below which that damping continues 
unabated. Further explorations of nonlinear Landau damping in the XY model are presented by Barre and 
Yamaguchi (2009). 



C. Smoothing the top-hat profile 

The XY model is convenient for exploring the differences between the dynamics underlying trapping and 
Hopf scaling. In this section we explore these differences in detail by considering the instabilities of a family 
of equilibrium profiles in which we smooth out the discontinuities of the top-hat profile of (32), As discussed 
in Sec. III.D, the top-hat profile (analogous to a broken-line shear flow or "water-bag" Vlasov equilibria) has 
a bifurcation to instability at a = 1 with c = 0, much like the Gaussian profile. However, we argued earlier 
(Sec. IV. A), the instabilities of the top- hat profile are expected to be affected by nonlinearity at strengths 
characterized by the Hopf scaling and be subcritical. In other words, the transition to instability for the 
top-hat profile is quite different from that for a smooth equilibrium profile like the Gaussian. 

The family of equilibrium profiles we adopt are given by 



R = 



tann — — — — tanli 



T 



T 



(172) 



where T is a parameter controlling the narrowness of the hyperbolic tangent functions. When T is decreased 
from order one values, the steps in the profile at / = ±1 sharpen, and for T limit to the discontinuities 
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FIG. 15 Scaling data for the first maximum, A*, plotted against distance to the stability boundary, 
a — ac, for the smoothed top-hat profiles of (172) for various values of the smoothing parameter, 
T. Panel (a) shows data for T = 1/16, 1/8, 0.26, .32, 1/2 and 1. Panel (b) plots the same data 
logarithmically (omitting T = 1/16), and includes additional data for T = 0.25 and 0.28. The 
dashed lines indicate the prediction of the single-wave model, and the dotted curve in (a) is the 
subcritical weakly nonlinear solution expected for the top-hat. 

of the top-hat. Using these profiles, we again conduct suites of initial-value problems; as for the Gaussian 
examples above, we begin these computations with p{0, /, 0) = and kick the system into action by adding 
the forcing (171) to the potential, ^{6,t), with taken to be small (10~^), at least initially. 

Figure 15 displays scaling data for the smoothed top-hat profiles of (172). For T > 1/4, the first 
maximum, follows the trapping scaling predicted by the single- wave model sufficiently close to the 
stability boundary (see panel (b)). However, as the profile becomes sharper, the trapping prediction for 
A* increases exponentially quickly (with A» proportional to (i?Q,)~^ oc cosh^(l/T)), and the system 
is unable to remain in the trapping regime much beyond onset. Instead, the saturation level flattens off 
quickly beyond the trapping regime and A^ approaches a different scaling that we interpret to be more 
like that for a genuine top-hat profile (e.g. T = 0.28). Indeed, for the lowest values of T (e.g. T = 1/8), 
the trapping regime is no longer observed and remains finite as a ^ Uc, which is consistent with the 
subcritical bifurcation of the top hat profile. 

The detection of a finite-amplitude branch of solutions for a > Gc with low T suggests that an unstable 
subcritical branch is born in the bifurcation at a = Oc, and that this branch subsequently turns around at 
smaller values of a, That is, there must be a second, saddle-node bifurcation for a < which generates 
the larger-amplitude stable solutions that return to higher a and provide the end-states of the initial-value 
computations. Thus, for the top-hat profile, or its relatively sharp smooth relatives of (172) with T -C 1, 
the instability has an abrupt onset. 

For intermediate values of T {e.g. T = 1/4), the system first follows the trapping scaling and, therefore, 
initially shows a relatively smooth transition. However, the branch's amplitude subsequently increases 
quickly and becomes vertical. The branch then disappears (no doubt by turning back on itself), leaving the 
system to grow to higher amplitudes and achieve a first maximum closer to that of the corresponding top-hat 
solution. The amplitude data in Fig. 15(b) is thereby rendered discontinuous. In other words, although 
the profile with T = 1/4 initially suffers a second-order phase transition following the trapping scaling, the 
transition can appear to become first-order if the system is pushed too far beyond the single-wave regime. 

Because we begin the initial-value computations from a weakly perturbed equilibrium profile parameter- 
ized by a, the solutions do not converge to any finite amplitude solutions for a < Uc, but Landau damp back 
to the unpatterned equilibrium. Nevertheless, it is possible to initialize computations with stronger kicks 
in this parameter regime (i.e. perturb the system with larger values for Aq). The resulting initial- value 
problems then display a clear threshold behavior, as illustrated in Fig. 16: for lower values of ^o, the sys- 
tem executes regular, weakly decaying oscillations with an amplitude set by ^o, as seen for the cases with 
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FIG. 16 Panel (a) depicts A{t), for a = 0.96 < Oc ~ 0.986 and T = 1/8, for computatfons kicked 
with different amplitudes for the forcing (171) added to ^{9, t). The values of are as shown 

(s) (u) 

and the "stable" and "unstable" first maxima, Al ' and Ar, are also indicated. For = 0.5, 
the stars indicate the times of the three snapshots of p{6, T, t) shown in panel (b) to the right. 

^0 = 0.4, 0.45 and 0.49 in Fig. 16. Simultaneously, p{9,I,t) develops no persistent phase-space pattern. 

(s) 

Larger forcings, however, trigger a growth of the mode amplitude upto a level, .4* , that is otherwise inde- 
pendent of as seen for = 0.5 and higher in the figure. The system subsequently executes trapping 
oscillations and p{d,I,t) develops a cat's eye pattern much like the end states generated in the linearly 
unstable regime. 

The threshold forcing amplitude Aq corresponds to a first maximum, Ai""^ that we assume is characteristic 

(s) 

of the unstable subcritical nonlinear state. On the other hand, the nonlinear states identified by ,4* , are 
plausibly the stable solutions on the upper branch. By recording the two first maxima, ^1"^ and Ai'\ we 
are therefore able to continue the scaling data of Fig. 15 to a < Oc for the cases with sharper equilibrium 
profiles (T = 1/16 and 1/8). As expected, the data for ^1*'' continue the original results for the first maxima 
from a > flc into the linearly stable regime. Moreover, this stable branch eventually terminates at lower a 
by coUiding with the unstable branch characterized by Ai^\ As illustrated in Fig. 15(a), the latter converge 
towards the weakly nonlinear solution branch predicted for a true top-hat profile (as given by the weakly 
nonlinear theory of Sec. IV. A. 3 which, after a translation of |^| = ^—2a2 back into the original unsealed 
variables, furnishes yll"'' — \/2{l — a)). 

In summary, top-hat (or broken-line and water-bag) equilibrium profiles display a very different nonlinear 
dynamics at onset to that of smooth profiles. The differences between the two in linear theory is well 
appreciated: the linear eigenspectra of smooth, stable profiles are characterized purely by the continuous 
spectra, which leads to Landau damping (or, more generally, to "continuum damping" for shear fiows, for 
which algebraic decay also appears), whereas piecewise linear profiles support persistent discrete waves riding 
on the interfaces of the equilibrium. That the weakly nonlinear states are also very different is not, however, 
so well acknowledged. 



IX. VARIANTS AND DEVIANTS 

So far we have seen several examples of single-wave model dynamics, viz., the patterns that form when 
equilibrium gradients inside the critical layer are wound up to saturate an unstable mode or arrest Landau 
damping. We next mention a number of variations on this single-wave theme. 



A. Long-wave theories and subharmonic instability 

An important ingredient in our derivation of the single-wave model is the imposition of periodicity in 
x, which leads to a wave-like global mode with a critical layer. Different boundary conditions would imply 
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FIG. 17 Subharmonic instability in an unstable jet. The panels show the late evolution of the jet 
after an unstable mode with x-wavenumber three grows to create a symmetrically placed pair of 
cat's eye patterns. In contrast to Fig. 6, the lower overtones were allowed to be excited by round- 
off error during the earlier parts of the computation, resulting in a subharmonic, vortex pairing 
instability. 

a different global mode structure in x and furnish a "single-mode" theory near the onset of instability. 
Alternatively, in spatially extended systems entire bands of wavenumber become unstable, rather than a 
distinguished neutral mode. In asymptotic analysis, one accounts for a finite bandwidth by introducing 
multiple lengthscales, with the long length scale describing the spatial modulation of unstable patterns. 

For systems suffering a long-wave instability the situation is more straightforward, since then one retains 
only a single, long length scale. When the unstable system also suffers critical-level singularities, the global 
mode amplitude depends on both the slow time and the long length scale, and the resulting weakly nonlinear 
theory couples familiar nonlinear wave equations to the critical-layer dynamics. Following this route in shear 
flows, Balmforth and Young (1997) and Stewartson (1981) coupled critical layers to the nonlinear Schrodinger 
equation and Boussinesq equations, respectively. More recent developments are provided by Shagalov et al. 
(2009). 

When the unstable wavenumber is finite, an alternative to solving the full multiple scale theory is to 
reconstitute the two spatial variables back into a single one. A single- wave-like model is then recovered, but 
contains additional spatial derivatives capturing larger-scale spatial modulations of the wave patterns. This 
route has been followed by Shagalov et al. (2010). 

Models with spatial modulation incorporate a key physical effect: the subharmonic secondary instability 
of a chain of cat's eyes that can lead to pairings of the vortices or plasma holes (e.g. Flierl et al., 1987 and 
Feix et al., 1994, respectively). This instability is easily illustrated for the simulation of a jet pictured earlier 
in Fig. 6. In that simulation, the third Fourier mode in x is unstable and grows to create the cat's eye 
pattern displayed in the figure; the lower overtones initially had zero amplitude, and any round-off errors 
that might otherwise excite these modes were deliberately suppresed. Once those lower overtones are allowed 
to develop, on the other hand, we arrive at the evolution shown in Fig. 17, in which a pairing instability of 
the cat's eyes emerges to destroy the pattern. This dynamics is not captured by the single- wave description, 
and is one of its main limitations. 



B. Modified single-wave models and forced critical layers 

The instability in the single-wave model results from the equilibrium gradients within the critical layer; 
once these gradients are twisted up by the global mode, the instability saturates (unless it can be reactivated 
by viscous diffusion). In other situations, the global mode can be rendered unstable by a different mechanism, 
and the instability competes with Landau damping in a critical layer (Balmforth and Young, 2002; see also 
Warn and Gauthier, 1989). The weakly nonlinear theory describing this situation is much like the single- wave 
model (with 7 = -1-1), except that the amplitude equation for A(t) contains a linear growth term. However, 
the unstable mode is not then able to saturate: as the mode grows, it twists up the critical-layer distribution 
into cat's eyes, removing the stabilizing equilibrium gradients that damp the mode. Consequently, the mode 
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amplifies even faster, twisting up an ever wider critical layer until the trapping scaling is broken. 

The opposite physical situation arises when the mode is driven by unstable equilibrium gradients within 
the critical layer, but damped globally (such as in the two-species plasma model with disparate ion masses 
explored by Balmforth and Kerswell (2002)). Then, the amplitude equation contains a linear decay term, and 
saturation again arises when the destabilizing equilibrium gradients in the critical layer are twisted up into 
cat's eyes. However, the global damping is unaflFected by those critical layer rearrangements and continues 
to damp the mode beyond saturation, forcing a subsequent decline in mode amplitude. The iiltimatc decay 
of the mode appears to be halted by the migration of the cat's eye across the critical region, which allows 
the mode to access fresh equilibrium gradients (Balmforth and Kerswell, 2002). 

Another natural variation on the theme is to add more waves to the problem (multi-wave models). In 
fluid mechanics this has been done for coupled Rossby waves (Vanneste, 1996), wind-excited water waves 
(Alexakis et al., 2004; Reutov, 1980), over-reflectional instabilities in shallow- water and compressible shear 
flows (Balmforth, 1998), and even-odd mode interactions in jets (Leib and Goldstein, 1989). In plasma 
physics, additional waves have been added in rclativistic beam-plasma systems (Evstatiev et al., 2003, 2005) 
and in models of plasma turbulence (Gary et al., 1992; Doxas and Gary, 1997). Other problems that could 
be attacked in this way (but have not yet been) include beat-wave resonance in plasma theory (Grawford 
et al., 1986), vortices (Mitchell and Driscoll, 1994) and astrophysical disks (Mattor and Mitchell, 1996), ion 
acceleration in the central plasma sheet of the Earth's geotail (Padhye and Horton, 1999), and the instability 
of inviscid flow over compliant walls (Benjamin, 1963; Landahl, 1962). All these cases involve at least one 
wave with a critical layer, together with wave- wave interactions captured within global amplitude equations. 

The critical layer in the single-wave model is rearranged due to the growth of a global mode. In other 
problems, there is no such mode and the critical layer is directly forced by an external perturbation. Models 
describing these situations contain a critical-layer equation in which the potential or streamfunction, tp, 
is determined not by an evolution equation but directly from the integral critical-layer distribution and 
the global forcing (as in the degenerate form of the single-wave model; Sec. V.B.2). The classic example 
of such a problem in fluid mechanics is the forced Rossby wave critical layer where waves excited by a 
uneven boundary enter a shear flow and impinge on the critical level (Redekopp, 1977; Stewartson, 1978; 
Warn and Warn, 1978). However, the theory has also been applied to a differentially rotating accretion 
disk perturbed by a forming planet (Balmforth and Korycansky, 2001), which is key in modeling angular 
momentum exchange and the migration of extra-solar planets. A novelty of forced critical layers is that the 
global potential or streamfunction contains a richer wavenumber content than the single-wave model. As 
a result, the filamentation of the critical layer distribution can suffer secondary instability (Killworth and 
Mclntyre, 1985), with the filaments rolling up into smaller vortices and generating more filamentation. A 
cascade ensues, leading to what Haynes (1985) has coined "critical-layer turbidcncc". 

All of the critical-layer models described above have the common feature that inside the critical region, 
the differential advection in y is locally linear. Theories that consider a critical layer located at the tip of 
a fluid jet where the flow is locally parabolic (Balmforth et al, 2001b; Brunct and Haynes, 1995; Brunet 
and Warn, 1990) lead to a shearing term with quadratic form y^(x- Disturbances in the linear initial- 
value problem are then found to decay algebraically as well as exponentially. Thus Landau damping is 
supplemented by a weaker form of continuum damping near the shearless point of the velocity profile. 
Moreover, the phase-space patterns that are generated can show a richer topology (see del-Gastillo-Negrete 
and Morrison, 1993), taking the form of either stacked cat's eyes or other kinds of patterns (Balmforth et al., 
2001b). 



C. Singular neutral modes 

Another key feature of the single-wave model is that the distinguished neutral mode that we perturb 
off the stability boundary is smooth. However, the distinguished mode need not be regular. For example, 
although symmetrical jets like that illustrated in Fig. 6 have smooth neutral modes, when the profile is 
made asymmetrical it generally becomes impossible to satisfy the condition, U"{y^) = U{y^) — Cr = 0, at 
both critical levels simultaneously. The neutral modes on the stability boundary then become singular, and 
the inner expansion for the critical layer is needed in the leading order of the asymptotics. Other examples 
include instabilities of Rossby waves (Hickernell, 1984), stratified and compressible flow (Ghurilov, 1999; 
Goldstein and Leib, 1989), and two-species plasmas (Balmforth and Kerswell, 2002; Berk et al., 1999). 
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Because^ oim; has to rc!-ordcr the asymptotic expansion for singular distinguished modes, the trapping 
scaling is no longer the relevant measure of the mode amplitude. It turns out that nonlinearity first becomes 
important at the level e'l'^ (a conclusion also reached by Crawford and Jayaraman, 1996, 1999, using center- 
manifold methods), and the analysis, which must be taken to much higher order than in the single- wave 
case, eventually predicts a delay-differential equation for the mode amplitude (Balmforth and Kerswell, 2002; 
Berk et al, 1999; Churilov, 1999; Goldstein and Leib, 1989; Hickernell, 1984). Unfortunately, that equation 
predicts that nonlinearity is destabilizing: beyond the initial phase of linear growth, the nonlinear terms 
promote an explosive growth of the mode amplitude (unless there is a significant viscosity as in Berk et al. 
(1999); Churilov (1999); and Goldstein and Leib (1989)). As a result, the mode rapidly leaves the "singular 
scaling" and only a transient is described by the model. 



X. CONCLUDING REMARKS 

For a pattern theorist, the single-wave model could be regarded as unsatisfactory because it offers no 
reduction in the dimension of the problem. Yet despite the lack of dimensional reduction, the single-wave 
model has several virtues over the original equations. First, mathematically, the model equations offer some 
simplifications over the full problem, but nevertheless captures the important dynamical behavior, such as 
the creation of cat's eyes. Second, the model appears in a wide range of different problems, and therefore 
provides a universal description of pattern formation in a variety systems with continuous spectra. 

Given this universality, it is worthwhile emphasizing the key properties of the original problem that are 
required for the single- wave model to emerge as the weakly nonlinear description in the vicinity of onset. 
As we have outlined, the single-wave model requires periodic boundary conditions in the direction of the 
equilibrium streaming or flow, and a regular neutral mode lying on the stability boundary to open the 
expansion of asymptotic theory. However, by themselves these properties are not sufficient to lead to the 
single-wave model: the so-called "Kuramoto model" describes a population of coupled oscillators with fixed 
amplitude and varying phase and natural frequency that is coupled through a mean field (Kuramoto, 1975; 
Strogatz, 2000). The continuum limit of this model furnishes a partial differential system that possesses 
equilibrium solutions corresponding to populations with randomly distributed phases. The associated linear 
problem has a continuous spectrum, stable populations exhibit Landau damping, and instabilities arise 
when a distinguished neutral mode detaches from the continuum (Crawford and Davies, 1999; Strogatz 
et al., 1992). However, the nonlinear dynamics of the instability, and the phase-synchronized population 
patterns that form, are quite different from the cat's eyes of our fluid and plasma problems (Balmforth 
and Sassi, 2000) . The similarity in the linear problem is therefore misleading, and simply the presence of a 
distinguished mode bifurcating from the continuum is not sufficient to anticipate the single-wave model. 

In fact, a key difference between the Kuramoto model and our plasma and fluid problems is that the 
former does not possess a Hamiltonian structure. An underlying Hamiltonian structure undoubtedly impacts 
any weakly nonlinear description. Nevertheless, the complete classification of models for which the single- 
wave model is an appropriate normal form has not been performed, and an underlying Hamiltonian structure 
may be neither necessary nor sufficient. There is, however, a Hamiltonian interpretation of the bifurcations 
to instability of the single-wave model. For example, the normal from case can be viewed a a continuum 
version of the Hamiltonian-Hopf bifurcation, which generalizes the discrete Krein (Krein, 1950) transition 
that is known to occur in a variety of fluid and plasma systems (e.g., Cairns, 1979; Casti et al., 1998; Kueny 
and Morrison, 1995; MacKay and Saffman, 1986; Morrison and Kotschcnrcuther, 1990; Sturrock, 1958; 
Tassi and Morrison, 2011; Tassi et al., 2008). This will be the subject of a future publication (Hagstrom 
and Morrison, 2013). 

Some other open questions arc as follows. First, the single-wave model describes the saturation of an 
unstable mode that twists up the destabilizing equilibrium gradients in its critical layer. However, in the 
single-wave scaling, the critical layer is much thinner than the range of unstable gradients, and saturation 
occurs when the background equilibrium has been flattened over only a small fraction of that range. Why 
does the mode not continue to grow and feed off the unstable gradients outside the narrow critical region? 
The answer surely lies in the fact that the system is constrained by the need to conserve the infinite number 
of Casimir invariants of Sec. II. A. Also, when dissipation is added to the system, the mode begins to grow 
secularly beyond saturation. Thus, the saturated single- wave system looks very much like a negative energy 
state of a finite-dimensional Hamiltonian problem with dissipation-induced instability (Krechetnikov and 
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Marsden, 2007; Morrison, 1998; Morrison and Kotschenreuther, 1990). 

Second, as mentioned earlier, it is not clear what is the ultimate fate of the unstable mode. Beyond 
saturation, aperiodic trapping oscillations occur that attenuate somewhat with time. However, one cannot 
decide whether the oscillations persist, or if the system relaxes to a steady BGK inodc\ One; often observes 
the intermittent creation and subsidence of secondary structures, or macroparticlcs, which contribute to 
the aperiodic motion (Tennyson et ai, 1994), albeit possibly with a dominant frequency. If the trapping 
oscillations continue indefinitely, then the critical region is in a constant state of agitation, and prolonged 
chaotic mixing occurs across the cat's eyes (del-Castillo-Negrete, 2000). 
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